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I. INTRODUCTION 

The breaking of chiral symmetry is an essential feature of the vacuum structure of QCD. It is a non-perturbative 
phenomenon that can only be studied from first principles by means of lattice QCD simulations. Because of the 
Nielsen-Ninomya theorem, exact continuum type chiral symmetry on the lattice necessarily implies the presence of 
doublers which have to be eliminated. Two widely used methods to do so are staggered fermions, which break the 
flavor symmetry but preserve a “chiral” symmetry, and Wilson fermions which preserve the flavor symmetry but 
explicitly break the chiral symmetry. 

An important ingredient of the vacuum structure is the topology of the gauge held conhgurations and its relation 
to the zero modes of the Dirac operator according to the Atiyah-Singer index theorem. Strictly speaking there is no 
topology at non-zero lattice spacing. For the Wilson Dirac operator one can obtain the topological charge by the 
spectral how in the quark mass [T]. Other ways to dehne the topological charge on the lattice is through the index 
of the overlap operator [2]as well as through spectral projectors 014]. One can equivalently employ purely gluonic 
methods where one measures a lattice discretized version of the volume integral of the topological charge density. The 
most common methods involve the gradient how [^, cooling 0]5] and APE/HYP smearing [ini E]. Of course all 
these methods agree up to cut-off effects but it is only the spectral how and the index of the overlap operator that 
give an integer result. We refer the reader to [12] for a critical comparison of all methods. 

The order parameter for chiral symmetry breaking is the chiral condensate which, according to the Banks-Casher 
formula, is directly related to the spectral density of the smallest Dirac eigenvalues. It is therefore important to have 
a detailed understanding of discretization effects on the Dirac spectrum and in particular on the topological zero 
modes. This problem was studied by means of chiral Lagrangians for the Wilson Dirac operator for QCD with three 
or more colors in the fundamental representation min]. Using mean held theory it was shown that the gap of the 
Hermitian Wilson Dirac operator closes when entering the Aoki phase. The microscopic Dirac spectrum was evaluated 
in great detail by means of a supersymmetric extension of the chiral Lagrangian and exploiting its equivalence with 
chiral random matrix theory [13 US]. One of the main results is that the Gaussian broadening of the topological zero 
modes scales as \fV'a = a for small values of the lattice spacing a which was hrst observed by lattice simulations 
m In the thermodynamic limit the real eigenvalues develop a band with a width proportional to VW^a^ = a^. A 
second important result is that the first order scenario |18j can only occur in the presence of dynamical quarks while 
in the quenched case at non-zero lattice spacing we necessarily have a transition to the Aoki phase when approaching 
the chiral limit m- We have also shown that the low-energy constants can be determined in a simple way by the 
properties of the Dirac spectrum [5D] . 

In the present article we consider the spectrum of the Wilson Dirac operator for QCD with quarks in the fundamental 
representation of SU(2). This theory reaches the conformal window for a smaller number of flavors as compared to 
three color QCD EH 122] and is of interest for modeling physics beyond the standard model, see for example two 
interesting reviews by Rummukainen and Kuti [231124j . What distinguishes this theory from QCD with three or more 
colors is that SU(2) is pseudo-real, and so it is possible to construct a gauge held independent basis for which the 
Dirac operator becomes real. For QCD at non-zero chemical potential this gives rise to a fermion determinant that is 
real so that the theory can be simulated for an even number of pairwise degenerate quark flavors, and is frequently 
studied as a model for QCD at non-zero chemical potential [25ll28j . More physically, this implies that quarks and 
conjugate quarks are in the same havor multiplet resulting in a doubling of the havor symmetry group. We expect 
that the discretization effects for two colors are qualitatively the same as for QCD with three or more colors. The 
main difference is the repulsion between the eigenvalues of the Dirac operator which is linear for small spacings. 

We compute the microscopic spectral density of the Hermitian Wilson Dirac operator starting from a supersym¬ 
metric extension of the chiral Lagrangian for two color Wilson fermions [2S] in the epsilon regime. The calculation is 
much more complicated than for QCD with three colors, and it is not clear whether the density of the complex and 
real eigenvalues of the Wilson Dirac operator can be obtained analytically. Moreover we compute the distribution 
p^ of chirality over the real spectrum of the non-Hermitian Wilson-Dirac operator and the mass dependence of the 
level density ps at the origin. Both quantities provide lower and upper bounds for the level density preai of the real 
eigenvalues of the non-Hermitian Wilson Dirac operator. We also consider the chiral condensate in the theory with 
and without dynamical quarks. 

Before calculating the observables mentioned above we briefly recall the properties of the Wilson Dirac operator 
for two colors and the definitions of the spectral observables in Sec. |ll] In Sec. jlllj we discuss the chiral Lagrangian 
for fermionic and bosonic quarks. The supersymmetric partition function for the quenched theory is evaluated and 
discussed in Sec. EYi Thereby we consider the continuum limit as well as the limit of a very coarse lattice and the 
related thermodynamic limit. The exact results and the detailed discussions of the single spectral observables are 
summarized in Sec. [V] Concluding remarks are made in Sec. jVIj In the appendices we briefly recall some properties 
of Bessel functions, present the detailed computation of the quenched partition function, and briefly rederive some 
spectral observables of continuum QCD with two colors and the Gaussian orthogonal ensemble. 
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Some of the results that appear in this article were first presented at the conference LATTICE 2012 |30) . 


II. QCD WITH TWO COLORS 


In subsection IIA we briefly recall the properties of the four dimensional QCD Dirac operator for the two color 
theory and fermions in the fundamental representation. Additionally we explain how chiral perturbation theory 
and the corresponding symmetry breaking pattern result from these properties. The variables we consider will be 
introduced in subsection III B I 


A. Dirac Operator for QCD with Two Colors 

For four dimensional QCD with two colors, the massless Dirac operator D anti-commutes with 75 implying chiral 
symmetry. Moreover, D commutes with an anti-unitary operator, m 

[D,At 2 C]= 0 , ( 1 ) 

where K is the complex conjugation operator, C = 7274 the charge conjugation matrix, and T 2 the second Pauli 
matrix acting in color space. Because {KT 2 CY = 1 it can be shown that there exists a gauge independent basis 
for which the matrix elements of D are real. In the Dyson classification of random matrix ensembles this symmetry 
corresponds to /3 d = 1 - 

Let us consider Ni quarks with a diagonal mass matrix m = diag(mi ,... Then the fermionic part of the 

Euclidean Lagrangian is given by 

r - ( \ I'o', 

) V o ijVf 1 o TO y ^”2 y ’ ^ 

where 


dp — -f f Ap 


(3) 


is the covariant derivative, (74 = —il 2 and Ca are the two-dimensional unit matrix and the three Pauli matrices, 
respectively, acting in Dirac space. The transformation property of the vector held under complex conjugation, 
A* = —T 2 ApT 2 is the reason for the anti-unitary symmetry Q. It also implies = —(J 2 T 2 <J^d^a 2 T 2 . This 

symmetry can be used to rewrite the two terms of the Lagrangian ([^, |25j 


C — Cq -\- Crm 


as 


and 


din — 


r = - 

m 2 


Cr 2 T 2'01 


ri 

(T2T2'4>1 


0 TOpdp 
TOpdp 0 


Cr2T2 


0 TO 
— TO 0 


l-JVf 


( <J2T2i^2 
V '^■2 


ri 

a2T2lpl 


iCTada ® I2N, 


f '02 

V Cr 2 T 202 


01 * 

CT 2 T 201 


If 02 

2 V f^2T'202 


(T 2 T 2 


0 TO \ / 'lp2 

-fh oil CT2T202 


Equation ([^ implies that the havor symmetry group is enhanced to U(2Af) given by the transformation 


(4) 

(5) 

( 6 ) 


( 0 r 

V 0-27201 


-)> 1 (g) C/* 


ri 

0-27201 


f 02 

V Cr2T2lp2 


1 g) [/ 


02 \ 

0-27205 J 


(7) 


with U S U(2Af). The havor symmetry is spontaneously broken by the formation of a chiral condensate which is 
given by the mass derivative of the partition function. Since the havor symmetry group does not act on a 2 and T 2 
in Eq. the symmetry is broken to Sp(2A^f) [^resulting in the symmetry breaking pattern U(2A^f) —>■ Sp(2A^f) 


^ The unitary symplectic group Sp{2A^f) is the section of the unitary group and the symplectic group which is sometimes also denoted as 
USp(2A^f). It should not be confused with the symplectic group itself which is non-compact and is sometimes denoted by Sp(27Vf), too. 
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[32L I33j . After taking into account the axial anomaly from configurations with a non-trivial topological charge v ^ 0, 
the Goldstone manifold is given by SU(2Aff)/Sp(2Aff) with 2Nf — Nf — 1 Goldstone bosons. 

When introducing the lattice spacing a in the Dirac operator Dw = D — aA one has to choose one of many fermionic 
discretizations. We are interested in the lattice artefacts of the Wilson Dirac operator, where the additional term in 
the Lagrangian ([^ is given by [34] 


/Iw = — 0i 


a 

2 


riVfA 0 

V'z; 'v 0 A 



f ri 

V f^2T2'lpl 


V 

I A(J2T2 G 


0 IjVt 
— IvTf 0 


ri 

(T2T2'01 


( 8 ) 


v-i Y 
0-2T2'lp2 J 


A(T2T2 


0 IaTj. 

— IWf 0 


V'2 \ 

f72T2^2 J 


Note that the four-dimensional Laplacian A = satisfies the symmetries A* = T 2 AT 2 and [A,(T 2 ] = 0. Hence the 
Wilson term ([^ transforms in the same way as the mass term Q. 


B. Spectral Observables 


The central object of our studies is the partially quenched partition function 

^ ^ ^ /det(-aA-|-7^(i^-bml-bxo75) \ 


(9) 


where (.)jyj is the average of two-color QCD with Nf dynamical quarks. The masses of the Nf dynamical quarks 
are encoded in the matrix M which may have also non-diagonal elements. The axial mass Xf has a non-vanishing 
imaginary part ie to guarantee convergence of the integrals. 

The generating function 0 allows us to derive two Green functions depending on whether we differentiate with 
respect to m' or Xf. The differentiation in Xf yields the Green function corresponding to the Hermitian Wilson Dirac 
operator D 5 = 75 DVV, i.e. 


Gf,{M,m,X,a) = ,Xo,Xi, a) 

Its discontinuity yields the level density of D 5 , 


m—m ^ 

Xq—X\—X 


tr—i— 

^ \ + Xl/ 


P 5 {M, m, A, a) = — lirn ImG 5 (M, m, A -I- ie, a) = — ~ 


( 10 ) 


( 11 ) 


Aff 


where X\(rn) is an eigenvalue of D 5 -\- 111^5 with respect to the eigenvector |fc). 

The level density p 5 satisfies an inequality with respect to the level density Preai of the real eigenvalues of the 
non-Hermitian Wilson Dirac operator This can be seen by considering the level density (11) at A = 0, 


P 5 {M,m,X = 0,a) = ^ (^ 6 [Xl{rn)]\ . 

\ k / Nt 


( 12 ) 


The corresponding eigenvalue equation reads 

(D 5 -b m 75 )|A:) = A|(m)|fc) 


(Dw - A|(m) 75 )|fc) = -m\k). 


(13) 


Hence the value —m is a real eigenvalue of Dw when A®.(m) = 0 which we denote by A^. Moreover taking the scalar 
product with the bra vector {k\ in the first equality of Eq. (13) and differentiating with respect to the mass fh yields 


9mAfe(m)|A6(iy,)=0 = 9m(A:|(D5 -bTO75)|fc)|;^6(^)^o = ik\l5\k). 


(14) 


The derivative of the vectors (fc| and \k) drops out because of the eigenvalue equations (fc|(D 5 + fh^^) = 0 and 
(D 5 + m 75 )|fc) = 0. The expectation value (fc| 75 |fc) is the chirality of the vector \k). Thus the slope of the spectral 
flow with respect to the quark mass at its zeros gives the chiralities of the real modes [T] . 
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Expressing the eigenvalues \l{m) in terms of and plugging Eq. ([I4| into Eq. (12) we find the relation 


P5(M,to,A = 0,a) = :^/^(5 [A|(to)]\ = ^ ( 

\ fc / Ni 


d[Af + m] 




= ^(i: 


Ni 




d[X^ + m] 

l(fc|75|fc)l 


.(15) 


N{ 


Since the modulus of the chirality of a normalized vector, {k\k) = 1, is always less or equal to 1 we obtain the 
inequality 


„(M,,„,A = 0,a) = -/ £ 


> ^ ( X! +m]) = Pree^liM, TO, a), 


V 


(16) 


Nt 


UWgIE 


Nf 


where the right hand side is indeed the level density of the real eigenvalues of Dw- 

The derivative of the partially quenched partition function with respect to the mass rh' yields another Green function 


G'{M,m,X,a) = dm'Z^(M,m,m',xo,xi,a) 


1 


m—m ^ 

Xq—Xi—X 


^ \ Dyf + toI + A75 ) 


= i(tr 


(17) 


Setting X = ie in the limit e —>■ 0 we can take the real part and the imaginary part. The real part is equal to the 
mass-dependent chiral condensate 


S(M, TO, a) 


1 

V 



1 

Dw -I- ml 



(18) 


In the case that to equals one of the eigenvalues toi, ... , 7 x 1 of the mass matrix M, e.g. say to = toi, one does 
not need the partially quenched partition function ([^ but the result can be immediately derived from the partition 
function 


a) = det(Ilw + • 


Then the chiral condensate is given by 




The imaginary part of G' yields the distribution of chirality over the real eigenvalues 


Py(M, to, a) =-lim Im G'{M, m, ie, a) 

TT e-aO 


1 


lim ( tr 


1 


1 


2TriV e-ao 

-lim ( trys 

2TTiV E^O ^ 


Dw + ml — iej 5 Dw + ml -(- zeys 

1 1 


-I- TOys — iel D^ + mys -I- iel 

= ^ (Xl'^Afc(™)]l(fc|75|fc)| ) • 


Ni 


N{ 


\ k / 

Here we employ the same relation between A^ and A^ as for the level density ps such that we have 


PxiM,m,a) = ^ {J^S[X^ -b TO]sign((fc|y 5 |fc)) \ . 

/ N{ 


(19) 


( 20 ) 


( 21 ) 


( 22 ) 


This expression explains the name of this quantity since it is the average sign of chirality at a fixed real eigenvalue 
Af = —TO. Understanding in this way has two implications. First the distribution of chirality over the real 
eigenvalues is normalized with the index u, 


/ dmp^{m) = ^ sign[(A:|y 5 |fc)] = v. 

AWgR 


( 23 ) 


















6 


We recall that the index v \s & topological invariant such that also this normalization is well-defined for each single 
configuration. It counts the total number of spectral flow lines that start at —oo for m —>■ —oo and end at -foo for 
TO —>• - 1 - 00 . Note that p^ifh) does not have to be necessarily positive definite. 

The second implication of Eq. (22) is another inequality with the level density of the real eigenvalues of D^^, 


\PxiM,m,a)\ = — 


+TO]sign((A:|75|fc)) \ 
\ k / 


Nt 


1 


<y{Yl = Preal{M,m,a). (24) 


\aW6IE 


Nt 


Combining this inequality with Eq. (16) we have an upper and lower bound for the level density Preab 


Px{M,m,a) < \p^{M,m,a)\ < preediM,m,a) < p^{M,m,\ = 0, a). 


(25) 


This inequality can be integrated over to such that we find an estimate for the average number of additional real 
modes for the set of configurations with a fixed index v which is 


0 < iVadd = -^real “ \k'\ < 


1 

V 


E 


1 - (fc|75|fc) 

|(A:|75|A:)| 


Nt 


(26) 


These estimates encode the only information available about the level density Preai as long as there is no closed 
expression for this quantity. Such a closed expression was indeed derived for QCD with three colors in [TCI ITC] but 
seems to be much harder to obtain for the two color case. 


III. CHIRAL LAGRANGIAN AND PARTITION FUNCTION 

We consider the epsilon regime of lattice QCD. In this regime the quark mass to, the level spacing a, the momentum 
of the Goldstone bosons p, and the amplitude of the Goldstone bosons 11 scale with the volume V as follows [55H57| 

TO ~ ~ ~ (27) 

In this particular scaling limit the Goldstone bosons with momentum p ^ 0 decouple from those with zero momentum 
such that it is legitimate to consider only the effective action of the latter. The corresponding chiral Lagrangian 
is space-time independent. Hence the integral (sum) over the lattice volume V can be performed yielding a prefactor 
of the chiral Lagrangian. 

It is clear that the epsilon regime is not the most physical regime. However, the mass scale of the epsilon domain is 
the one of lowest eigenvalues which are most sensitive to lattice artefacts. Using the leading order epsilon domain chiral 
Lagrangian one can also determine the region where the artificial phase structures like the Aoki-phase m manifest 
themselves, see also subsection |IVD| Another important application is the measurement of physical parameters such 
as the chiral condensate, the pion decay constant as well as the electroweak effective couplings from lattice simulations 
in the epsilon regime [38114(1] . 


A. Ghiral Lagrangian for the Fermionic Quarks 


The chiral Lagrangian for QGD with two colors follows from the requirements that the effective partition function 
should have the same flavor transformation properties as the QGD partition function in the full theory. Similar to 
three color QGD the partition function in the epsilon domain is given by an integral over the Goldstone manifold 
U(2A^f)/Sp(2iVf). A matrix U € U(2Af)/Sp(2Af) is an antisymmetric unitary matrix parametrized hy U = VIV^ 
with I defined by 


/ = 


0 Iatj 
— Iat, 0 


(28) 


The chiral partition function for a fixed topological charge ly is given by 


Zf^(M,a) = 


l\J{2Nt)/Spi2Nt) 


d/r(U)det'^/Vexp ^tr {MIU + {MI)'<U-^) - ahi {{lUf + [lU)-^) 


(29) 
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Since the integrand is a class function, the integration can be extended to U(-/Vf) with the integration measure given 
by the Haar measure. Without restriction of generality we can assume that > 0 is non-negative. To keep the 
notation simple we have introduced the dimensionless quantities, 


M = VT.M and = VW^a^ 


(30) 


Thus the low energy constants E (chiral condensate) and Wg are absorbed in these new quantities. Moreover we wish 
to recall that the low energy constant Wg can be a priori positive or negative such that the effective lattice spacing 
a may also appear as an imaginary number even though the physical lattice spacing a is always real and positive. 
However, as is the case for fundamental Wilson quarks for QCD with three or more colors, we find that Wg > 0 in 
the case of two colors. This is discussed in subsections |III B| and |V A] The general mass matrix MI is antisymmetric, 
but we restrict ourselves to a mass matrix such that M = (m + ixQ)l 2 Nt proportional to the identity matrix. Here we 
have to distinguish between the real quark mass m and the real axial mass xq. The latter is introduced as a source 
for generating particular observables, see subsection |H B| 

The chiral Lagrangian in Eq. (291 is not the most general lowest order Lagrangian in the epsilon regime. We can 
add the terms VWg(Pti'^I{U — U^) and VWra^ti^I{U + U~^), which are of the same order as the other terms. As 
in the case of three color QCD these two terms have the same symmetry properties as the Wg term, and they have 
to be included a priori. However, based on large Nc arguments [41], we expect that the strongest effect results from 
the Wg term. Moreover the Wg and Wr term can be easily introduced by a convolution with a Gaussian in the quark 
mass m and in the axial mass xq, respectively, see |5D|. Therefore we stick to the partition function Z^^{M,a). 

The coset of antisymmetric unitary matrices is equivalent to the circular symplectic ensemble CSE(2Aff) which is 
the ensemble of self-dual unitary matrices. Matrices in this ensemble satisfy the relation C/^ = —lUI and are related 
to antisymmetric unitary matrices via the scaling U —)■ lU. This set of matrices was first studied by Dyson [44j . 
Thus the eigenvalues of U are pure phases and are Kramers degenerate. In terms of the CSE(2A^f) the chiral partition 
function (29) reads 


Zl^^{M,a)= [ d/r(f7)deWVexp hr {MU + M^U-^) - ahr + U-'^) 

JCSE{2Ni) [2 


(31) 


which may be more familiar to some readers. The measure dii{U) on the Goldstone manifold CSE(2A^f) is the 
normalized Haar measure induced from the Haar measure on U(2A^f) and can be extended to an integration over 
'U{2Nf) as in case of the representation (29). The measure can also be expressed in terms of the flat measure dU 
which is the product of differentials of all independent matrix entries, i.e. the measure is |45j 


dfi{U) oc 


dU 


det^'-i/V 


(32) 


It is invariant under the group action U —>■ IV^’^IUV for any V G U(2Nf) because the Jacobian of this transformation 
is d(IV'^IUV) = det^^^~^VdU. When diagonalizing U = 5'diag(e*‘^E..., ) 0 125'^ with S G Sp(2A^f) the flat 

measure reads 


Nf 

dU(x (e*‘^'= - U‘^‘fdn{S) 

l<k<l<N{ fc=l 


(33) 


with dfi{S) the normalized Haar measure of the coset Sp(2A^f)/Sp^'(2). 

To interpret the term in Eq. (31) as a diffusive process, see [46] for three color QGD, it is quite convenient to 
express also this term as a mass term convoluted with a Gaussian, 


exp [-a^tr {U^ + = exp [-a^tr {U - U-^ f - 


(34) 


= (7^(6 / dcrexp 

Jse\i{2Ni) 




with the constant = 2'^f/^(8o^7r) ^f(2JVt i)/ 2 ^ matrix a is self-dual and Hermitian a = = —la"''I and the 

set of those matrices is denoted by Self(2A^f). Thus the partition function (311 simplifies to 


zO^*(M,a) = / daf dpiU)deU^^U exp 

Jse\l(2Nr) JcSE(2Nr) 

Z^^{M + ia,a = 0), 


{{M + ia)U + {M + ia)^U ^)--^trtT^ 
z Iba^ 


= CNtC / dcrexp 

JseU{2Nt) 


1 

16a^ 


trcT 


(35) 
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which is expressed in terms of the partition function at a = 0, see |42L I46j . This expression has some interesting 
relations to the microscopic partition function of the chiral Gaussian orthogonal random matrix ensemble (chGOE) 
which is up to a factor the same as the partition function Z^^{M + ia,a = 0), 

^chGOEd-^ + ia][M^ — id]) = det~^^^(M^ — + ia,a = 0), (36) 

and of the Gaussian orthogonal random matrix ensemble (GOE) of finite matrix dimension v x ly, 


z: 


GOE 




f dadet^^^{M^ — ia) exp 

Jse\i(2N{) 


1 


16a2 


trcr 


(37) 


partition functions, i.e. 


In subsection IV B we will argue that in the continuum limit the partition function (351 factorizes into these two 

(38) 




Thus the v eigenvalues of the finite dimensional GOE lie on the scale a and shrink to a single Dirac delta function at 
the origin showing that they are the former v zero modes of continuum QGD with two colors. 

Since the matrix a is self-dual its eigenvalues are Kramers degenerate. The flat measure da can be rephrased in 
terms of these eigenvalues via a diagonalization a = SsS"^ = S'diag(si,..., satJ G with S G Sp(2iVf) yielding 

Nt 

da oc n {Sk - si)'^dfi{S) dsk- (39) 


When additionally assuming that the mass matrix M = (m + ixQ)'k 2 N{ is proportional to the identity matrix the 
diagonalizing matrix S drops out in Eq. (351. Note that S can be absorbed by U because U is Haar distributed. Then 
the partition function (35) reduces to 


/ I'oo r c2 1 \ 

Z^^{M,a) j dsk{m + ixQ + iskY exyp j {sk - siY + ixo]l 2 Ni + is). (40) 


Therefore the chiral partition function at finite lattice spacing can be understood as a convolution of each single quark 
mass with a Gaussian. Since these Gaussians are not completely statistical independent we can understand such a 
convolution as collective fluctuations as we have already seen those for three color QGD, cf. [min]- 


B. Chiral Lagrangian for Bosonic Quarks 


To find the bosonic chiral Lagrangian, we start from the observation that it is possible to find a gauge independent 
basis in which the massive Hermitian Wilson Dirac operator of Vb bosonic flavors, = ysDw = 757 /id/i G liVb + 
j 5 Z>fh+lZ>x — 075 A 0 lAfbj is real and symmetric because of the anti-unitary symmetry. We need to generate its 
inverse and, hence, its matrix Green function such that the axial mass x must have a non-vanishing imaginary part 
which we assume to have the signature L — diag(-|-lp, —IjVb-p) with 0 < p < iVb positive signs. Then the inverse 
square root of the determinant of can be represented as an integral over real bosonic fields {(jyf, ^|’), 


Y det Ds 


= C 


J dcj)id(j)2 Gxp i J ( 


d'^a 


</>i 

4>2 


1 ® Lx -I-1 1 




Lm — aA i 


L 




iL 


1 ® LX — 1 0 Lm -I- oA 0 L 


<(>i ^ 

h J 

(41) 


with C a normalization constant. The multiplication of the matrix iL guarantees the existence of the integral. 

To eliminate the square root we have to double the number of flavors IVb —>■ 2fVb. The Dirac term of the Lagrangian 
then reads 


Lq = - i{a^dY*) ® L(j)2, (42) 

with L = L ® I 2 where the two-dimensional unit matrix reflects the doubling of the number of flavors. This term is 
invariant under 


-> 1 G U(j)i, 4>2 ^ 1® [L{Lf ^)'^]4>2 


(43) 
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with U € Gl(2iVb,M) because the bosonic fields are real. With maximum breaking of chiral symmetry [351 [53], the 
set of transformations that leaves the “bosonic” chiral condensate, 


E oc 



1® L 0 
0 -1®L 


</'2 


^ 0 , 


(44) 


invariant reduces to the set of matrices satisfying U^'^LU = L and U G Gl(2iVb,K). This set is known as the 
non-compact orthogonal group 0(2p, 27Vb ~ 2p). The Lorentz group 0(3,1) is a particular example of this kind of 
non-compact groups. 

The pattern of spontaneous symmetry breaking is given by Gl(2A^b,K) —t 0(2p, 2A^b — 2p) and the resulting 
Goldstone manifold Sym^(2p, 2A^b — 2p) = Gl(2fVb,K)/0(2p, 2A^b — 2p) is the set of real matrices U € K2Arbx2Afb 

the symmetry [/^ = LUL and 2p positive real eigenvalues and 2A^b ~ 2p negative real eigenvalues. Such matrices 
can be parametrized hy U = LWW'^ with W G Gl(2A^b7R)- Another parametrization is via a diagonalization 
U = Odiag(e*b ■ ■ • I with O G 0(2p,2Nt — 2p) and Sj G K. The corresponding invariant 

measure dfj,{U) is given by 

dll 

dp{U) (X -—’ with dG (X (e®'=-|-e®') je'*''— e'*'| d/r(0) de'*'" 

^ l<k<l<2p l<k<2p 2p+l<k<l<2N^ k=l 

2p+l<l<2N^ 

(45) 

and dp{0) the invariant measure on 0(2p, 27Vb ~ 2p). Hence dpiU) is indeed invariant under the group action 
U —?> LV"^LUV with V G Gl(2iVb,K) because the fiat measure (product of the differentials of all independent matrix 
entries) transforms as d(LV^L[/V) = ^VdU. 

Let us first discuss the bosonic partition function for p = 0 or p = 2fVb. In these cases the imaginary part of x 
has only positive or negative entries, i.e. L -G with L = ±1. Then the non-compact orthogonal group reduces 

to the compact 0(2fVb) and the coset is equal to the set of all real symmetric positive definite matrices Sym_|_(2fVb). 
There are no integrability issues, even for a = 0, since the convergence of the integrals is assured by the imaginary 
part of X. For a ^ 0, the convergence of the integral is already guaranteed by the term if p = 0, 2iVb. 

For p = 0, 2A^b the chiral partition function is given by 


f dp{U)det exp 

(^trTO([/ - U-^) + ^tr x(17 -b t/-^) - ahi (H^ -b U-^) 

/Sym^( 2 Afb) 

2 2 


(46) 


The sign L = ±1 does not drop out although the imaginary part of x can be dropped for finite a. Because the 
eigenvalues of U are positive definite, the transformation U -G —U is not allowed. However we can absorb the sign 
in front of the mass by changing U -G U~^ which yields a change of the topological charge v —)■ —v. There is some 
freedom in the choice of the signs and phase factors of the terms in the action, but the sign of the term is fixed 
by the physics of the problem. To have convergent integrals for ^ 0 we necessarily need that the sign of this term 
is negative. We also require that the bosonic and fermionic terms in the action can be combined into a supertrace. 
This is achieved by rotating the fermionic and bosonic variables by a factor i and including an overall minus sign in 
the bosonic action which gives the representation of Eq. (46), cf. the fermionic partition function (311. 


For p ^ 0, 2N}, the integral over 0(2p, 2Nt — 2p) is non-compact such that the measure (451 cannot be normalized 
even by group invariant potentials. We always have to include an infinitesimal increment in x, even for ^ 0, in 
order to get a convergent integral. The bosonic partition function is then given by 


Zy{m,X,a) = 


f dp{U)det exp 

l-tr mL{U—U ^) + Irtr {xL + ie){U + U ^) — a^tT{U^ + U ^) 

' Sym^ (2p,2Arb-2p) 

2 2 

(47) 


C. Supersymmetric Partition Function and Quenched Theory 

Having discussed the bosonic and fermionic partition function we are now ready to formulate the supersymmetric 
partition function that generates the microscopic Dirac spectrum. Especially we focus on the average of the Green 
function which can be traced back to derivatives of the partially quenched partition function ([^. Hence we put 
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iVf — >■ iVf + 1 in Eq. (31 1 and set A^b = 1 in Eq. (|46|). Then the combined chiral partition function reads 
Z,{M,X,a) = C-^ I 


dfi(U)Sdet'^^'^U exp 


i:(2Aff+2|2) 


with 


M = 


MO 0 
0 ml2 0 
0 0 m'l2 


-^Str {MU - M'fU-^) - ^Str !([/ + U-^) + o^Str {U^ + IJ-^) 


(48) 

(49) 


X = 


0 0 0 
0 a;ol2 0 
0 0 a;il2 


and 
C = 


f dfj,{U)Sdet'^^^U exp 

-^Str {MU - M'fU-^) - ^StrX(17 -h 17"^) -h a^Str {U^ + G'^) 

'^{2Nt+2\2) 

2 2 J 


m—m 

Xq^Xx 


(50) 

The mass matrix of the dynamical quarks again satisfies M = —IM"^I. The normalization constant guarantees that 
the partition function is equal to 1 if m = to' and Xg = Xi because the determinants in Eq. cancel. A supermatrix 
U in the Goldstone manifold E(2iVf + 2|2) = U(2Xf + 2|2)/UOSp(2Xf + 2|2) has the block structure 


(U, -Ig^\ 

\ g 


(51) 


where U{ G CSE(2Xf + 2), Ug G Sym^(2), and g is a 2 x (2Xf + 2) matrix whose matrix entries are independent 
Grassmann (anti-commuting) variables. The whole matrix fulfills the symmetry 


U^ = 


Uf -g^ 

gi 


— M2N{+2 0 
0 I 2 


Gf -Ig^ 

g Ub 


4l27Vf+2 0 

0 I 2 


— 4l2Wf+2 0 

0 I2 


u 


4l2Wf+2 0 

0 I2 


(52) 

We recall that “T” acts on a supermatrix as the supertransposition which has a slightly different action on the 
off-diagonal blocks (the sign of the upper right block changes). Moreover, we choose the following notation of the 
supertrace and the superdeterminant 


StrC/ = tr Gf — tr [/b, Sdet 17 = 


detjUt + iIg^Uy^ g) 
det 17b 


(53) 


This choice differs with respect to some other works SZlllH] where the sign of the supertrace is reversed and the 
superdeterminant is the inverse of the definition here. The invariant measure dg,{U) is given by 


dg{U) = 


dUfdUhdg 

Sdet^f-i/2c/’ 


(54) 


where dUf, dUg, and dg are the flat measures meaning the product of the differentials of all independent matrix 
entries. It naturally appears in superbosonization |49H51j which can be used to directly map random matrix theory 
to a finite dimensional version of the chiral partition function [52] . 

Eor the quenched theory where N{ = 0, the partition function drastically simplifies 


Z^{M,X,a) = C-^ j dUfdUbdgSdet^^+'^'^/^Uexx, 

S(2|2) 


-yStrM(C/ - t/-i) - yStrX(t/ -h C/-^) -h a^Str (17^ -f C/-^) 


where we already wrote the measure in terms of the flat one. The matrices in this expression are 


M = 


TOI 2 0 
0 to'12 


X = 


^ f a;ol2 0 
I 0 a;il2 


, U = diag(l2,0) 


0 a* 

0 


ci-V - 


a a 


a —/? 

0 


diag(l2,0’ 


(55) 


(56) 


/3 /3* 0 e® 


with O G 0(2). The orthogonal matrix O drops out and yields a constant. Thus the remaining integration is over 
a phase with Lp G [—tt, tt], two positive variables and with Si,S 2 G K, and four Grassmann variables 
a, a*, /3, /3*. The measure becomes 

dUfdUbdg -G dU^ |e®i - e^^\de^Me‘‘^dada*dl3dl3* = to'‘^+®i+®=^ le®^ - e^^\dipdsids2dada*dpd/3*. (57) 
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We recall that the integration over Grassmann variables is defined as f da = 0 and f ada = 1 and similar for the 
other Grassmann variables. Integrations of higher orders are not needed to be defined since Grassmann variables are 
nilpotent because of their anti-commuting nature. 

As for the partition function of fermionic quarks (31) we can linearize the term in by introducing an auxiliary 
supermatrix a, 


Z^{M,X,a)=C-^ [ da j d/r(C/)Sdet'"/^C/exp 

Jt.(2\2) Jy.{2\2) 


E(2|2) 


X exp 


-yStrX(C/-hC/-^) 


/S(2|2) 


1 

16a^ 
d[a] exp 


StrCT^ — ^Str (M — ct)(17 — [/ 


1 

Iba^ 


Str a^ 


(58) 


The set S(2|2) consists of (2|2) x (2|2) supermatrices which can be parametrized as follows 


a = diag(l 2 , 0 ) 


iu 0 rj* X* 
0 iu —r] —X 
r] T]* vi 0 
X X* 0 ^2 


diag(l2,0'^) 


(59) 


with O S 0(2), u,vi,V 2 G R and ??,?7*,X;X* independent Grassmann variables. Additionally, the supermatrix a 
fulfills the symmetry ( [5^ . With the help of the partition function of the chiral Gaussian orthogonal ensemble, 

(60) 


^ihGOE([^ +X-a][M-X-a]) = Sdet -’'^‘^{M-X - a)Z^{M - ct, X, a = 0), 


we can rewrite the quenched partition function as a convolution of the ^chGOE "with a Gaussian 


Z,{M,X,a) = 


dtrexp 


16a2 


Str (T^ 


c{a = 0 ) r 
C{a ^ 0) yg( 2 | 2 ) 

z''X’oo^{[M + X - a][M - X - a]), 


Sdet'^/^(M-X-cr) 
1 


/S(2|2) 


d[a] exp 


16a^ 


Str ( 


(61) 


This representation shows the diffusive character of the effect of a finite lattice spacing a on the spectrum of the Dirac 
operator, see [S3] for the diffusive approach in supersymmetric spaces. Thus it exhibits a similar effect as already 
found for lattice QGD with three colors of the Wilson Dirac operator HOlllSl. Interestingly, when omitting the term 
■^chGOE integral (61) we obtain the finite dimensional quenched partition function of an orthogonal Gaussian 

random matrix ensemble. 


z(^) 

^GOE 


M-X 

4a 


/S(2|2) 


dcrexp 


1 


16a^ 


Str 


Sdet 


•'^^{M-X-a)/ [ d[a] 
/ JE(2|2) 


exp 


1 


16a^ 


Strtj"^ 


(62) 


We recall that the real eigenvalues of this ensemble scale with the lattice spacing a which can be also seen in this 
expression. The normalization constant C is independent of a because the partition function should be normalized for 
any choice of M = TOI 4 and X = XI 4 , cf. Eq. (|^ where we started from. Thus the ratio is C{a = 0)/C{a ^ 0) = 1. 
Now we have all ingredients to calculate and discuss the quenched partition function. 


IV. QUENCHED PARTITION EUNCTION 


In subsection |IV A| we present the result for the quenched partition function at finite lattice spacing, quark mass 
and axial mass. We reduce this result to a compact integral which represents the fermionic valence quark and a 
non-compact two-fold integral which reflects the nature of the bosonic valence quark. The continuum limit o —t 0 
is derived in subsection IV Bj The limit of a very coarse lattice |a| 3> 1 and the related thermodynamic limit are 
discussed in subsections |IV C| and |IV D| respectively. 


A. Discussion of the Partition Function at finite a 


The evaluation of the partition function ( |55| ) is straightforward but tedious. 
Grassmann variables 


Since we have only four different 
the Grassmann integrals can be evaluated by a brute force expansion. This is worked out in 
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appendix B 1 The final result for the generating function is given by 


ZAM.X.a) = \r^ 


dsi 


ds -2 


sinh 


Si - S2 


^iy{2iip Si S2)/2 g^p 


X exp [—2iLaio cost/j + iLxi{coshsi + coshs 2 ) + 4a^ cos2i^ — 2a^(cosh2si + cosh2s2)] 

X [(40^(6“^*“^ + + iL{m — Xo)e“*'^ + iL{m' — a;i)e“®^ — — 1) 


x(4a"(e 


2/'p-2i¥> 


c, —2S2 


_l_ gj<P+s 2 _|_ g iLim — X{))e + iL{m'— xi)e — z/— 1) 


+40^(36-2*^ + + e-^®i + + 26-*“^-®! + 26-*“^-®") 

+iL{2{m — xo)e~^‘^ + (m' — a;i)(e“®^ + e“®^)) — i/ — l] 

= 16a'^(<i)ji_45'ii_0 + ^iySiy+4^0 + ^iy+2Siy-2,0 ~ ^iz-2iS^+2M + 4<i)ji_25'i/+2,2 + 24>ii_iS'i,_i4 
+24’,i_3iS';i+i 4 + 8<i)ji_|_i5',y+i_3 — 6$,i+iS'ji+i 1 + 2 $ji_iS',i_|_3 1 + 4$ii5',i_ 2 ~ 2$jiS'ji o) 

+4a^((2z2 — l)$,y_2iS'i/,o ~ (2 z 2 + l)4’ji5';i+2,o + 4z/$,yS'ji_|_2,2 + 2(i2 + l)<i>,i_|_iS',i_i_i + 2(z2 — l)<i>,y_iS'ji3_i^i) 

—8a^(m — a;o)(4>ii_3S'jy_o + 2$ji_iS'iy+2,2 — 4>,y_i5'i,+2,o + + <i)i,_2S'i/+i,i) 

—8a^{m' — a:i)(4>,i_2S';i+i^i + $jiS'ji _|_3 i + 2 <i)ii+iS'ji ^2 ~ 4>;i+iS'ii^o + 4>;i_iS',y+2,o) 

+{m - xoy^v- 2 Siyfi + 2{m - xo){m' - - 2z/(m - xo)^u-iS^fi 

+{m' - xi)‘^^^S„+ 2 ,o - 2i2{m' - + (z 2 + (63) 

Each term factorizes into a non-compact two-dimensional integral S'^^a and a compact one-dimensional integral 
which are defined by 


4*- 


/ TT ^ 

exp [2Lm sin ip — 2iLxo cos p + 4a^ cos 2p\ 

-■K 27’" 


= 


dp 

2n^ 


{m? — 

(m -I- xq)^ 


OO 


1 — — 00 


exp [2771 sin p — 2ixo cos p + 4a^ cos 2p\ 
h (40^) 1^+21 (^2^ 


xq — m 
Xo m 


— xt 


(64) 


and 


S^^aim' ,xi,a) = 


'ay 


dsi 


dS2 


sinh 


Si - S2 


o“7‘(si+S2)/2 


cosh^ 


Si - S2 


X exp [7Lrz7'(sinhsi -I- sinhs 2 ) + iTa;i(coshsi -I- coshs 2 ) — 2a^(cosh2si -|- cosh2s2)] 

X exp [7m'(sinhsi -|- sinhs 2 ) + 7a;i(coshsi -|- coshs 2 ) — 2o^(cosh2si -|- cosh2s2)] 

/ OO nOO 

ds / dye~^^y°‘ exp [2iLm'y sinh s -1- 2iLxiy coshs — 4a^(2y^ — 1) cosh 2s] , (65) 

-OO J 1 


where we have performed the substitution y = cosh(si — S2)/2 and s = (si -|-S2)/2 in the second line of Eq. (65). The 
function /; is the modified Bessel function of the first kind. 

The partition function is correctly normalized which can be readily checked either by setting oY = m = m! = 0 and 
xq = Xi = iLe with e —)■ c» or by setting m = m' = Xq = xi = 0 and oY ^ oo. In the latter case where we take 
a —>■ OO we have to distinguish between even and odd v because of the following asymptotics, 


<i)^(m = 0, Xq = “>■ 00 , a = 0) 


Sfi^ai'fn = 0, ccq = ie —>■ oo, a = 0) 


{-aye^’^ 

^ 2v^ ^ 

yY(iLYe-^^ 


$^(777 = 0, Xq = 0, a —)■ 00 ) 


•S';!, 0,(777 = 0, Xq = 0, a —)■ 00 ) 


{ayd 


li„2a-‘ 


mod2(i^),05 


ysY a 
v^(-7L)'^e-2“ 


( 66 ) 


2£3/2 ’ '--7 2ll/2a3 

The partition function should be also unity when m = m!, xq = xi, and a finite. However this general situation is 
quite non-trivial to verify. It goes back to Cauchy-like integrals of superfunctions invariant under supersymmetric 
groups [ISJ [MH55] which were first derived in a general form by Wegner . Indeed the integrand in Eq. (55) is 
invariant under the supergroup UOSp(2|2) when m = m! and xq = xi. 

We emphasize that the representation of the compact integral (64) as a sum is highly convergent for 777 , Xq, oY small 
enough since the Bessel function behaves as 7/ oc 1/d for large order. This cannot be said about the non-compact 
twofold integral (651 which causes some trouble if the lattice spacing a is too large or too small or when the quark mass 
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m! or the axial mass xi are too large. Then we have large cancellations which challenge the numerical integration. 
Hence it is advantageous to improve this integral. Luckily we are either interested in its value for xi = 0 or in its 
imaginary part. In appendix |B 3| we calculate both expressions and find 


= 0, a) = 


[sign(A^)]^ 

— 1 )! 


dy 


lyCH+l 

dt — , exp [—4a^(2j/^ — 1)] 


'o v/ 22/2 - 1 


(67) 


and 


t * 


1 


|t' + t|l^lexp 


8a2(2j/2 _ 1) 


{t' + t — m'Y 


Kfj, + t|) — {f —>■ —t} 


t'=0 


Su,a{m', xi, a) = — lim Im Su,a{fn', Xi + ie, a) 

TT e—>-0 


( 68 ) 



= -(-sign^)^ 

“ (ImI - l)!fc! 


dy 


dt 


y 


,a—|/i| + l 


8a2(2y2_l) 


{t — sign{fi)m'Y — 4a^(2j/^ — 1 ) 


fc=0 


TT /.oo a+i 

X \ 2yY\ - c) ( 


8a2(2y2 _ l) 


{t — m'Y — 4a^(2?/^ — 1) 


Xi +1 
xi — t 


0(x^ — Y). 


The Heaviside step function 0(xf — Y) restricts the integral over t to a compact interval in the second term while in 
the first term the integral can be simply evaluated by a Taylor expansion in t and has to be skipped for y = 0. For 
fj, = 0 the term 1/t and the derivative have to be omitted and the minus sign in front of {t —> —t} becomes a plus 
sign in Eq. (67). The functions and Kf^ are the Bessel function of the first kind and the modified Bessel function 
of the second kind. 


B. Continuum Limit (|a| <C 1) 

To understand what happens when taking the continuum limit |a| —> 0 it is useful to consider a random matrix 
theory which is equivalent to the same chiral partition function as two color QCD with Wilson fermions. This random 
matrix ensemble consists of random matrices of the form 


Dyj 


aA W \ 
-W^ aB ) 


(69) 


with A = and B = real symmetric matrices of size nxn and {n+i')x (n + i/), respectively, and W an n x {n + i/) 
real rectangular matrix. The random matrix ( |69[ ) belongs to the symmetry class 19QC^ in the classification scheme 
of non-Hermitian matrices by Magnea m- When we choose the Gaussian 


P{Dw) = exp 


1 

16 


(tr + tr B^) - ^tr WW^ 
2n 


(70) 


as the distribution of Bw, we do not have to unfold its spectrum, i.e. in the limit of large matrix dimension n —>■ oo 
the smallest eigenvalues of Bw and B 5 = 75 BW (with 75 = diag(l„, —l„+jy) in its matrix form) around the origin 
are of order 0(1). We have chosen the Gaussian distribution for simplicity. Due to universality other probability 
distributions give rise to the same chiral Lagrangian. 

The random matrix (69) can be brought into the form 


B(y — 


Of 0 
0 Of 


B C 0 \ 
B' f ] 

orb) 


Oi 0 \ 

0 O 2 J’ 


(71) 


where A, B' , f and b are of 0{a), Oi G 0{n), and O 2 G 0{n + v). At vanishing lattice spacing a the matrix B(^ has 
v zero modes and n pairs of imaginary eigenvalues zLiy. This ensemble is known as the chiral Gaussian orthogonal 
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ensemble misiisi]. For small lattice spacing |a| ^ 1 we can apply first order perturbation theory. To lowest order 
in a, the secular equation factorizes as 


det(Z)w - Xl 2 n+u) = det 


7 -Al„ 

c Y 

.V 

-Al„ J 


det(6 — Alj,). 


(72) 


Therefore to leading order the effect of a non-zero lattice spacing is the broadening of the zero modes to the spectral 
density of the ensemble the matrix b is drawn from. By the central limit theorem, the distributions of the matrix 
elements of b are Gaussian even if the matrix elements of B do not have a Gaussian distribution. Hence the distribution 
of 6 is a dimensional GOE. The variance of b is equal to the one of B so that the distribution of b is given by 


p{b) oc exp 


r N ,1 

■ 

-Tvytrl)^ 

4a2 

= exp 


V 


32a^Wf 


■trb^ 


where we used the relation between the dimensionless lattice spacing a and the physical one when identifying N 
i.e. 


a^VWs = Na^/8. 


For v = 1 the distribution is a simple Gaussian 


p{b) = exp 


V 


32a^W8 


(73) 




(74) 


(75) 


In applying the central limit theorem, we have assumed that the average of the matrix elements of O 2 BO 2 vanishes. 
This is not the case if the average (Aki), (Bki) ~ aSki which indeed is the situation in lattice simulations, e.g. see 
[62j . However, such term is exactly of the form of a mass term and can be eliminated by a redefinition of the mass 
justifying that the central limit theorem gives a centered Gaussian (731. 


Indeed the splitting of the spectrum into to one of chGOE and one of GOE can be also seen at the level of chiral 
perturbation theory which also applies to the epsilon domain of QCD. 


Gonsidering the quenched partition function (55) we have shown that it can be rewritten as the convolution of the 


partition function of continuum QCD without zero modes and a Gaussian function of the supermatrix tr, see Eq. (611. 


Rescaling cr —>■ au we can drop the a dependence in the partition function ^chGOE([-^ + X — aa][M — X — acr]) since 
the function remains finite at a = 0. Hence the quenched partition function factorizes for |a| ^ I, 


Z^{M,X,a 


yiA 

~ ^chGOE 


(m^ - X^^ 


^GOE 


M - X 
4a 


(76) 


Both partition functions and their resulting spectral observables are summarized in appendix This factorization 
is a also shown in Fig. The comparison of the approximation (76) with Monte Carlo simulations of the random 
matrix model (69) and the full analytical result (63) confirms that this factorization applies quite well for |a| < O.I. 


Only for small indices, namely = 0,1,2, the deviations are persistent. The reason is the non-analytic nature of the 
spectrum of the continuum QCD Dirac operator at the spectral edge A = ±m ^ 0. The first derivative of the level 
density of the Hermitian Dirac operator diverges at the edge for z/ = 0,1, 2. Thus the limit a —)■ 0 cannot be 
interchanged with the integral over the eigenvalues when averaging over the spectrum implying that the continuum 
limit is not uniform in these cases. If the quark mass vanishes, only the case z/ = 0 is non-uniform. We return to this 
point in subsection |V B| 
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— v=0, Pi 

V=0, PchGOE+PGOE 

— V=l,p5 

V=l, PchGOE+PGOE 

— V=2, Pi 

V=2, PchGOE+PGOE 

— V=3, Pi 

V-3, PchGOE+PGOE 

— V=4, Pi 

V-4, PchGOE+PGOE 









FIG. 1. Comparison of some spectral observables of the Wilson Dirac operator with lattice spacing VM4D3 = |o| = 0.0625 
and of the spectrum consisting of the sum of continuum QCD without zero modes and a v x v dimensional GOE for indices 
{v = 0,1, 2, 3, 4). Figure a) emphasizes that the distribution of chirality py{m) over the real eigenvalues (solid curves), see 
Eq. p^ , and the level density ps (m, A = 0) of the Hermitian Wilson Dirac operator at the origin as a function of the quark 
mass (dashed curves), see Eq. (B17l, become the level density of a GOE. In subsection VE we argue that then also the level 
density of the real eigenvalues of the non-Hermitian Wilson Dirac operator Dw shares this distribution. The splitting of the 
spectrum of D 5 into the superposition of a chGOE and a GOE is shown in figure b) where the quark mass is set to m = 0. 
The level density ps (solid curves) perfectly agrees with the sum pchGOB + Pgoe (dashed curves) when the eigenvalue A stays 
away from the origin while the deviations remain close to the origin. Especially the case of vanishing index u = 0 shows 
obvious deviations which are discussed in subsection |V B| These deviations are particularly large for small index v — 0,1,2 and 
become more prominent when the quark mass m becomes non-zero, see figure c) for m = 2 (black vertical line). The spectral 
discontinuities of the continuum Dirac operator are hardly suppressed by the GOE level density pgoe shifted by the quark 
mass (dashed curves). However the level densities P 5 (m = 2, A) (solid curves) are smooth at any finite values of the lattice 
spacing. Only for larger indices, here iz = 3, 4, the agreement with a splitting into the two sub-spectra is almost striking, cf. 
figure d). The bad, see figure e), or good, see figure f), agreement carries over to the mass dependent chiral condensate E(m), 
see Eq. (18 1 . Also in these two plots the solid curves are the exact results (B18 1 at finite lattice spacing while the dashed curves 
are the approximation of the chiral condensate as the sum of the chiral condensate of continuum QGD without the contribution 
from the z ero modes, SchGOE, and of the “chiral condensate” resulting from a finite dimensional GOE, Egoe, see Eqs. (C7| 
and (|G30|). 
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C. Limit to a very Coarse Lattice (|a| » 1) 


We also consider the limit of a coarse lattice |a| ^ 1 to get a full understanding of lattice artefacts in the spectrum 
of the Wilson Dirac operator. Especially this limit shows the change of scales in the spectrum when increasing the 
lattice spacing. It may help to estimate the strength of the lattice artefacts. 

In the limit of a very coarse lattice the quark mass m and the axial mass x are an order smaller than . Then the 
compact integral (641 is dominated by the term 4a^cos2^. Expanding around the two saddle points ip = 5(p/a and 
Lp = TT + dpi a the function <I>^ becomes 


$^(TO,a;o,a) 


|a|>i 


2 'Ka 

i-^LY 


dSpe 


— SSip^ 


exp 


2Lm 


d(p — 2iLxQ 


+ (-l)^exp 


2Lm 


Sp + 2iLxQ 


exp 


4a^ + 


160^ 


^-2iLxo _|_ 


(77) 


A similar expansion can be done for the non-compact double integral (651. However this time we have only one saddle 
point namely Si /2 = 


5',x,Q(m',a:i,a) ' « exp[ 2 zLa:i - iY] 


8 Y 

[iLYY^ 


dSsi / dSs 2 l^si — ^S 2 | exp 


iLm! 


(i5si -I- 5s2) - 4((5si -I- Ssl) 


exp 


2iLxi — AY — 


m 

16Y 


(78) 


Since m,x Y we can omit all terms in the sum (631 which do not come with an Y factor. Using the approxima¬ 
tions (771 and (78 1 the quenched partition function becomes 

2iL(xi — xo) + 


Z^{M,X,a) « exp 


16a2 


(79) 


The second term in Eq. (771 depending on the sign (—1)^ cancels in the sum. Therefore the dependence on the 
topological charge is completely lost. 

We wish to emphasize that this scaling is only valid around the origin and applies for the eigenvalues which are 
of order 0{a), cf. Figs. A more physical scale for the level density of the real eigenvalues and the “mesoscopic” 
spectrum around the origin of the Hermitian Wilson Dirac operator is discussed in the next subsection. 


D. Thermodynamic Limit 


In this section we discuss the thermodynamic limit sometimes also referred to as the mean field limit where the quark 
mass m, the axial mass A, and the squared lattice spacing are of the same order, and satisfy m = YdVm, A = EUA, Y = 
WgVY 3> I. This limit is best performed in an eigenvalues representation of the supersymmetric integral (551. 
Choosing V = —iLU as a new variable, the exponent in the partition function (55) is given by 


exp 


I, 


I 


--StrM(U -h y-Y - -Str A(U - U"^) - YSti {V^ + U"^) 


The corresponding saddle point equation reads 


p{v) = jiv- y-Y + ^(^ + - 2a2(^^ - y-Y = 0. 


The solution of this equation is computed in appendix For A = 0 it is 


y = 



-iu 1 


-iY)\ 

signml 4 . 


I 4 , \m\ < SY, 

ItoI > SY. 


while for A 7 ^ 0 it takes the quite complicated form 

y = e '^'^^^'^14 with p = arccos 


m cosh d + sinh d 


8 Y cosh 2d 


(80) 


(81) 


(82) 


(83) 
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FIG. 2. In plot a) we show three different cases for the thermodynamic limit of the level density ps. The Hermitian Wilson 
Dirac operator D 5 exhibits a spectral gap around the origin when the quark mass \m\ is larger than 8 a^. This gap closes at 
m = 8 a^ and the system enters the Aoki phase. Hence the value of the level density ps at the origin is an order parameter 
for the Aoki phase. In plots b), c) and d) we compare the thermodynamic limit (black curves) of the chiral condensate, the 
distribution of the chiralities over the real eigenvalues, and the mass dependence of the level density ps at the origin A = 0, 
respectively, with the behavior at finite lattice spacing \/VWsa = a = 1. Although the finite a result has still large deviations 
from the thermodynamic limit, the phase transition building up at m = ±8a^ is clearly visible. Also the dependence of the 
observables on the index u has almost disappeared. 



I 


^2-A2 


1/3 


(8c 


212 


- I 


mX 

( 8 ^ 


' m?\^ 


I 


m? — A 2 
“( 8 ^ 


- I 


1/3 

. (84) 


The level density of the Hermitian Wilson Dirac operator can be calculated from these saddlepoint solutions via 


P5(m, A,a) = 


—ImStrH + —Strdiag(0, 0,1, l)Im (H — y ) 

ZTT 2 


(85) 


m—m 
X\—X2 — X 


The first term is sub-leading such that it can be omitted in the discussion The density ps is plotted in Fig. [^a) and 
d). It exhibits two different scenarios, either has a spectral gap at the origin or the gap is closed. The critical 
points in the case of a gap can be read off the solution (84) which are 


A = ±„|l-| t 

m 


2\ 2/3\ 3/2 


( 86 ) 


Instead of the involved derivation of these critical points via the saddle point solution presented in appendix one 
can also simply find them by a substitution W = signTO(y — V~^)I2 in Eq. (80) such that we have to minimize the 
function 


q{W) = |m|Vl4-biy2 + AIF - 


(87) 
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At the critical points the first two derivatives of qiW) vanish, i.e. 


q'{W) = 


mW 


VI 4 + 


+ A - 8a^W = 0 and q"{W) = 


m 


(l4 + Vh2)3/2 


- 8 a^ = 0 . 


( 88 ) 


This yields W = — 1 implying the critical point ( | 86 [ ) from the first derivative. Hence the spectrum of 

ZI 5 has a gap for m > 8 a^ in the interval 


A G 


8 a' 

m 


2 \ 2/3' 


3/2 


-HI- — ,H 1 - — 


8 a^ 


2\ 2/3' 


3/2' 


(89) 


For TO < 8 a^ and A = 0 the gap is closed and the system is in the Aoki phase. 

The distribution of chirality over the real eigenvalues and the mass dependent chiral condensate are given by the 
saddlepoints solution V via the relations 


Px(m,a) = 
E(to, a) = 


V. „ .. 1 . 


— ^ImStrH — Astrdiag(0,0,1, l)Im(H + y ^) 

27r 2 


,, 1 , 


— HR-sStrH — vStr diag(0,0,1, l)Re (y + y ^) 

27r 2 


m—m 

X1—X2—O 


m—m 

X1—X2—O 


(90) 


Since the axial mass is set to zero we only need to consider the solution (82). The first term is again sub-leading for 
I](to, a) while it is leading for a)- This different behaviour hints to a separation of scales which indeed happens, 

see Fig. 

We discuss the thermodynamic limit of the spectral observables in more detail in the next section. 


V. SPECTRAL STATISTICS OF THE WILSON DIRAC OPERATOR OF TWO COLOR QCD 


In this section we discuss the spectral observables in more detail and derive their analytical expressions. In par¬ 
ticular, we summarize the continuum and thermodynamic limit. We study the fo llowin g spectral observables: The 
unquenched partition function and its chiral condensate in subsection V A the level density ps of the 
Hermitian Wilson Dirac operator in subsection VB the chiral condensate Yi{m) of the quenched theory in sub¬ 
section [WCl the distribution of chirality Py over the real eigenvalues of the non-Hermitian Wilson Dirac operator Dw 
in subsection VD[ and the level density of the real eigenvalues and the number of the additional real eigenvalues in 
subsection IV e1 


A. Partition Function of Dynamical Quarks 


For Nf = \ there is no spontaneous symmetry breaking and the QCD partition function for QCD with two colors 
is the same as for QCD with three or more colors. Indeed for Nf = 1 the partition function (291 can be written as 


^(TO,a)= [ exppTOCOS p — 4a^ cos2ip], 

27r 


(91) 


which coincides with the one-flavor partition function for /? = 2 , see na. _ 

For two or more flavors the partition function is the Nf{2N{ — l) dimensional integral (29). For simplicity we assume 
that all quark masses are equal to to. Then the partition function 


Z^'{m,a) = f 

Jc 


CSE(2Art) 


dp(17)det‘"^^C/exp ^tT{U + U ^) — a^ti {U'^ + U ^) 


L 2 


(92) 


can be calculated by diagonalizing U and applying de Bruijn’s integration theorem [63] ■ Then we end up with a 
Pfaffian (essentially an exact square root of a determinant of an antisymmetric matrix) |64| 

Zl^^{m,a) ocPi[A] with Au = {k - = 0,a)A < kj < 2Nf. (93) 

where <i)g is the integral given in Eq. (64) evaluated of xq = 0, 

4)° (to, a:o = 0, a) = (—f)'^ f exp [2TOsinp-|-4a^ cos2p] . 

J-TT 27r 


(94) 
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FIG. 3. Chiral condensate for two dynamical quarks for index = 0,1, 2, 3 and lattice spacing x/VWga = a — 0, 0.125,0.25, 0.5. 
The continuum limit (black curves) diverges as vim for m —>■ 0 due to the zero modes. For larger index v, the peaks reminiscent 
of this singularity are more persistent at finite lattice spacing. The reason is the localization of real eigenvalues around the 
origin which are smoothed out only very slowly. 


The chiral condensate is easily calculated by taking the derivative of the free energy In (m, a) with respect to the 
quark mass resulting in 


S^f(m,a) = ^dmlnZ/^‘(m,a) = ^tr A ^dmA. 


(95) 


The derivative of A is given by 


dmAki = ik-l) [<^°_ 2 Nt+k+i im,xo = 0,a) + ^l_ 2 Nt+k+i -2 im,xo = 0,a)]. 


(96) 


The chiral condensate is normalized to the asymptotics 
for particular indices zz = 0,1, 2,3. 


lim 
2 —>- + 00 


i) = Nf. In Fig. 


we illustrate the behavior 


As already mentioned in subsection III A the unquenched partition function (92) factorizes in the partition function 


of continuum QCD without the zero modes and in the partition function oi a. vxv dimensional GOE in the continuum 
limit |a| 1. The spectrum of the finite dimensional GOE is of the order a which shrinks to the origin in the continuum 

limit. This creates a singularity of the chiral condensate at o = 0 for configurations with v 0, cf. Eig.[^ 

In the thermodynamic limit the partition function is either dominated by the saddlepoint Uq = signm for \m\ > 8a^ 
or by the saddlepoints satisfying Uq + for \m\ < 8a^. Then the resulting chiral condensate is given by 


W)jVf = ^E^f(TO,a) 
Nf 


m 

signm, 


|m| < 8a^, 
\m\ > 8a^. 


(97) 


This result does not depend on Nf and we will show in subsection |V G| that it is also valid for the quenched theory. 
The thermodynamic limit is already well approximated for a lattice spacing y/VWsa = a « 2, see Fig. 

Finally, let us show that the sign of Ws has to be positive. We start from the observation that the partition 
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FIG. 4. Comparison of the chiral condensate with two dynamical quarks for index u — 0,1, 2, 3,4 with the thermodynamic limit 
(black curve). Although around the origin the thermodynamic limit is a good approximation for lattice spacing \/VWsa = a — 1, 
we have still large deviations in the phase transition region around the qnark mass m = 8a^. At \/VWsa = a = 2 the 
thermodynamic limit is almost approached. Nonetheless the u dependence is barely visible already at smaller lattice spacings. 


function (92) has to be the same as the average over the eigenvalues over + 77175 which are real, 


Z^'{rn,a) = / ]J(A®(777)) 


Nt 


(98) 


£>6 


Hence the partition function has to be positive definite for Nf even. The sign (—1)'^'^* in front of the average results 
from multiplying 75 with the non-Hermitian Wilson Dirac operator D^. As is the case for /3 = 2 the partition function 
satisfies the general relation (see Eq. (29)) 


Z^^{m = 0 , 0 ^) = = 0 ,- 0 ^). 


(99) 


Therefore, the partition function for Nf G 4No + 2 cannot be positive definite for both values of the sign of Wg. This 
can also be seen from the explicit expression for the two-flavor partition function, which, at m = 0 is given by the 
two dimensional integral 

Z^^=-^{Tn = 0,a) = \ [ -e*‘^=|^e*"(‘^i+‘^^^exp[-4a2cos2(^i-40^0032(^2]. (100) 

2 277 J_^ 2tt 


All odd powers of the phases and vanish since the a dependent term has double the frequency. Therefore, 
for odd V, only two terms remain in the expansion of the Vandermonde determinant. They can be combined as 

Z^^=^{m = 0, o) = -4/(,,_i)/2(-4a2)/(,+i)/2(-4a2) = 4/(,_i)/2(4o2)/(,+i)/2(4o2). (101) 


This term is always positive definite if a is real and thus Wg is positive while it is negative when 
Hence a negative value of Wg contradicts to the positivity of the partition function. 

The unquenched partition function (101) is shown in Fig. for v = 1 and Nf = 2,4. Indeed the 
for Nf = A does not depend on the sign of at m = 0 (see Eq. but the figure shows that the 

is not positive definite for a negative value of Wg. 


VWgo^ = < 0. 


partition function 
partition function 


B. Level Density of D5 


The level density pg given by 


pg{m,X,a) = da;^lmZ^{M,X,a) 


M—ml4 

A = (A-|-i£)l4-i.Al4 


imm ediat ely results fr om c ombining Eqs. (10), (11), (63), and (BIO). We obtain a quite complicated expression (see 
Eq. ( |B17 ) of appendix B4| which can be evaluated numerically. In Fig. we compare this result with Monte Carlo 
simulations of the random matrix theory (69). 





























21 




FIG. 5. The mass dependence of the two flavor (left) and the four flavor partition function (right) for VWga^ = aF' = ±0.01 
and 1 / = 1. In the case where the low energy constant is negative we have always regions where the chiral partition function is 
negative although it is positive when considering the full theory. 


When taking the continuum limit a —>■ 0 the peak around the quark mass m shrinks to a single Dirac delta function 
with area v. At non-zero a, the peaks are no longer degenerate and broaden to a width ~ a. Due to the weak level 
repulsion which is linear for real matrices we do not observe separate peaks for v > 1. 


In the limit a —?► 0 the density of the cluster of v peaks around m combine into the level density pgoej see Eq. (C33l, 
oi a. V X V dimensional Gaussian, cf. Figs. 0b). c) and d). In this limit we expect that pa can be approximated by 

l“l<i |A| ( /7x-x\__ Ts fi/i (m-\ 






m- 


rPchGOE (\/A2 - m2) 0(A2 - m^) ± 


4a 


with Pqoe tbe level densities (C34), (C35), and (C33) depending on v and 


PchGOE(^) “ J dyJu{y^ ± 2|a;| (Jj,(2|a;|) Jy_i(2|a;|) Jy+i(2|a;|)) 


( 102 ) 


(103) 


the microscopic level density of continuum QCD and, thus, of chiral GOE without zero modes [65]. This approximation 
is shown in Figs. 0b), c) and d). 


As we already pointed out in subsection IV B the continuum limit via the approximation (102) is not uniform at 
X = m for small values of v. In particular, as can be seen by a brief computation (see appendix for m = A = 0 
and 1 / = 0 we never reach the correct value at the origin in the continuum limit, see Fig. [Tj Thus we obtain the limit 


lim p 5 (m = 0, A = 0, a) = —/= ^ 1 = lim P 5 (m = 0, A, a = 0). 

a->0 y'2 A->0 


(104) 


The non-commutativity of the two limits has also been confirmed by Monte Carlo simulations, see Fig. Therefore 
the continuum limit is not uniform in m so that one has to be careful with quantities which depend essentially on 
eigenvalues close to the origin. 

As was discussed in Eq. (15), the mass dependence of P 5 (m, A = 0, a) is given by the density of the real eigenvalues 


of Dw weighted by the absolute value of the inverse chirality of the states. In the continuum limit the density of the 
real eigenvalues is well approximated by the level density pooE , see Eq. (C33), of exactly the same finite dimensional 
GOE which also describes the former zero modes, cf. Fig. |^with height that is of order 1/a and a width that is of 
order a. Since the chirality of the states with real eigenvalues is |(fc| 75 |fc)| « 1, in the continuum limit, we also find 
that the mass dependence of P 5 (m, A = 0, a) is given by a GOE in this limit. 

In the thermodynamic limit the level density ps takes the form 


P5 


m A 

8a2 


cosh 




m A 


sin 




m A 


0 


A^ - - (8a^)^/^)' 


(105) 


where p and F) are given by Eqs. (83) and (84). The Heaviside theta function implies a spectral gap if |m| > 8a2. If 
the gap is closed the system is in the Aoki phase. The order parameter is the pion condensate which is proportional 
to dx log.^|x=o and, hence, to p^ijn, A = 0, a). In the thermodynamic limit this quantity becomes a semi-circle 




m A 


= 0 




(106) 
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FIG. 6. Comparison of the analytical result (B17l for the level density ps (solid curves) with Monte Carlo simulations (MC, 
symbols). We have simulated for two different values for the valence quark mass, m = 0 (blue curves) and VT,fh = m = 2 
(red curves). The index has been chosen v = 1 (left plots) and = 2 (right plots) and the lattice spacing is ^yVWsa — a = 
0.125, 0.25,0.5. For each plot we have generated 10® matrices distributed according to Eq. (701. The prominent peaks at m = 0 
and m = 2 result from the broadened former zero modes. They will become Dirac delta functions in the continuum limit. 


(see Fig. [^d). It immediately shows that the phase transition at |m| = 8a^ is of second order as in the case of three 
color QCD with Wilson fermions. Furthermore we can say that the height of the limit (106) is of order 0(1) and its 
integral is of order O(a^). These two pieces of information become important in the discussion of the real eigenvalues 
in subsection IV El 


C. The Chiral Condensate 

The analytical result for the chiral condensate 

T,{m,a) = dm'^eZ^{M,X,a) 


M—m]l4 
X—iel-4^0 


(107) 










































23 



FIG. 7. Continuum limit of the level density ps of the Hermitian Wilson Dirac operator. The index and the quark mass are set 
to zero, iz = m = 0. The continuum limit (black solid curve) is not uniformly approached by the analytical result dBlTt . We 
have confirmed our analytical result by Monte Carlo simulations of the random matrix model (691 and (70l. The ensemble of 
2.5 X 10® matrices gives an accuracy of about 1%. The dotted black curve is the limit for a very coarse lattice where ps = 2/7r 
is a constant. 




FIG. 8. The chiral condensate for parameters used for the simulation in Figs. |^and|^ We underline that the chiral condensate 
is for the quenched theory while the one shown in Fig. is for the theory with dynamical quarks. Interestingly, the chiral 
condensate in the quenched theory approaches the asymptotic value 2 from above while it is approached from below for the 
theory with dynamical quarks. 


for the quenched theory is explicitly shown in Eq. (B18) of appendix B4 Its behavior is shown in Fig. 

As in the theory with dynamical quarks we have a l/m singularity for m —>■ 0 due to the zero modes at non-zero 
index zz 7 ^ 0 in the continuum theory. These singularities are washed out by the finite dimensional GOE according to 
the approximation 


E(m, a) ' « -f (m), 

where the continuum limit of the chiral condensate without the contribution of the zero modes, S 


(108) 


chGOE("i): is given 


in Eq. (C9). The continuum result is in agreement with earlier work by Damgaard et al. | 66 j . Indeed the contribution 

works quite well for small but finite lattice spacing and \iy\ >2, cf. Figs, [^e) and f). The 


of see Eq. (|C30) 


large deviations for smaller indices result from the non-uniform continuum limit of the smallest eigenvalues. The 
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former zero modes help to push the spectrum away from the origin. This is the reason why the approximation by the 
GOE and continuum QCD without zero modes is very accurate for larger indices. 

The thermodynamic limit of the chiral condensate (B18) can be easily taken via the Eqs. (90) and (82) yielding 


(V'V')Aff=0 



m 

4^’ 

2 signml 4 , 


\m\ < 8a^ 
\m\ > 8a^. 


(109) 


This limit is already visible for a lattice spacing ^JVW%a = a 
obtained for dynamical quarks. 


1. Moreover, it is exactly the same result (97) as 


D. The Distribution of Chirality over the Real Eigenvalues 


Also the distribution of the chiralities over the real eigenvalues 

Py (to, a) = -9mlm Z^{M, X, a) 

TT 


M—ml4 

A=(A+ie)l4-i-Al4 


( 110 ) 


can be written as a sum of products of a compact integral and a non-compact twofold integral ( see Eq. (B19) of 
B4). Comparisons with Monte Carlo simulations of the random matrix model (69) are shown in Fig. 

In the limit a —>■ 0 the chirality distribution, py(A) becomes the level density pgoe of a z/-dimensional GOE which, 


for 1 / = 1, is a Gaussian, see Eqs. (C33), (C34) and (C35). We compared this behavior of with pgoe and p^{X = 0) 
in Fig. j^a) at a lattice spacing ^jvW^'a = a = 0.0625. The agreement is almost perfect. An interesting distinction 
from the case of three color QCD, where the distribution of the chiralities over the Dirac spectrum in the continuum 
limit also agrees with a finite dimensional Gaussian random matrix model, is that the peaks corresponding to the zero 
modes are only barely visible humps. As in the case of three color QCD the number of these humps is equal to ly. 
However they merge with the main peak at the origin due to the very weak level repulsion. This becomes apparent 
when comparing Figs. [^a) andj^with Fig. 1 of m- _ 

The thermodynamic limit of Py follows from the first equality of Eq. (90). The supermatrix V + V ^ is real at 
the saddle point (l82|) such that the second term in Eq. (90) vanishes. Hence the first term is the leading one. The 


derivative in the mass yields an inverse square root behavior, 

, . V 0(8a^ — ItoI) 

Py(TO) « -- 


TT \/ (8a^)^ — m? 


(Ill) 


Therefore, py of two color QCD exhibits the same square root singularities at the boundary of the support as in the 
case of three color QCD, see [5D]. We observe this asymptotic behavior in Monte Carlo simulations where it is starting 
to build up, cf. Fig. and in the numerical evaluation for even larger lattice spacings a = 1 shown in Fig. |^c). 


E. Level Density of the Real Modes of Dw 


Since the analytical derivation of preai for /3 = 2 was a tour de force [HI [20] , it is not surprising that we did not 
succeed to find an analytical result for the density of the real eigenvalues of the non-Hermitian Wilson Dirac operator 
Dw for two color QCD with (3 = 1. Nonetheless we can extract some information via the inequality (25) and thus 
via the distributions P 5 (to, A = 0) and Py(A). 

For small level spacing |a| <C 1 the distributions P 5 (to, A = 0) and Py(A) are well approximated by the density 
PGOE of the finite dimensional GOE, see Fig. [Tja). Therefore also the level density of the real eigenvalues has to 
be preai ~ PGOE for |a| <C 1. This is indeed the case already for the lattice spacing y/VWsd = a « 0.1, see 
F ig. Surprisingly the distribution of the chiralities over the Dirac spectrum agrees quite well with preai even for 
y/VW^a = a ~ 0.25 when the index is ^ > 2. Since the integral of the difference Preai — Px equal to the average 
number of additional real modes N^dd we deduce that this number is highly suppressed for configurations with a 
larger index ly at small lattice spacing. For the Wilson Dirac operator of three color QCD the average number A"add is 
of order see [Hj- The reason is that the probability of finding an additional real eigenvalue is proportional 

to 




(112) 
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FIG. 9. Comparison of the analytical results (B17l (red solid curves) and (B19l (blue solid curves) of the mass dependence 
of the level density P 5 (A = 0) at the origin and the distribution of chirality over the real eigenvalues, respectively, with Monte 
Carlo simulations of the random matrix model (691 (MC, symbols). We have generated an ensemble of 10® matrices both for 
index v = 1 (left plots) and index v = 2 (right plots) at lattice spacings V^VlFga = a = 0.125,0.25, 0.5. We have also calculated 
the distribution of the real eigenvalues preai of the non-Hermitian Wilson Dirac operator Dw (green symbols). As correctly 
predicted by the inequality (251 the level density preai lies always in between and P 5 (A = 0). Interestingly, for small lattice 
spacings, the level density is better approximated by Px than by P 5 (A = 0) for small lattice spacings. Moreover, we notice a 
significant difference between the integral of the three distributions for large lattice spacings while their support remains the 


where the real eigenvalues Xk are of order a. In the present two color case we expect that fVadd is of order C)(a'^+^) 
because of the smaller level repulsion. Our expectation is confirmed by Monte-Carlo simulations of the random matrix 


model (69), see Fig. 10 


In the thermodynamic limit the height of Preai has to lie between the square root singularity of on the scale 
0{a~‘^) and the semi-circle of P 5 (A = 0) on the scale 0{1). From the simulations shown in Fig. we notice some 
kind of flattening of Preai and a departure from for larger values of a. Moreover, the average number of additional 
real modes extracted from Monte Carlo simulations, see Fig. suggest a ly independent behavior where IVadd «: a. 
Since also the distribution of the real eigenvalues lives on the support [—8a^,8a^], we conjecture that preai becomes a 
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FIG. 10. Log-log plot of the dependence of the average number of additional real modes on the lattice spacing a. The Monte 
Carlo simulations (symbols) of the random matrix model (691 for index v = 0,1,2 were htted by the behavior oc a" for small 
lattice spacing (colored lines) and oc a for large lattice spacing (black line), separately. We generated random matrices for 
various matrix dimensions and ensemble sizes to keep the statistical and systematic error below 2%. The transition region, 
VFWsH = a = 0.25 — 2, for the scaling behavior is excluded from the fits. 


constant flat plateau on this interval with a height of order 0{a in the thermodynamic limit. Interestingly this is 
exactly same behavior as in three color QCD, cf. [201 • 


VI. CONCLUSIONS 

We have computed and analyzed the microscopic spectrum of the Wilson Dirac operator for QCD with two colors 
and the quarks in the fundamental representation. The discretization effects are very similar to the three color case. 
Especially in the thermodynamic limit, lUWga^l = |ap oc |UEto| = |m| oc |UEA| 1, they are essentially the same. 
When the quark mass is large enough, |to| > 8 a^, the microscopic spectrum of the Hermitian Wilson Dirac operator 
D 5 develops a symmetric gap around the origin of a width — ( 8 a^)^/^]^^^. The system enters the Aoki phase 

when the gap closes, or in terms of the eigenvalues of the non-Hermitian Dirac operator Dw, when the mass \m\ hits 
the strip of eigenvalues. In this phase, parity is broken spontaneously with a non-vanishing pion condensate that is 
proportional to the eigenvalue density of at the origin. In the thermodynamic limit, the support of the distribution 
of chirality over the real eigenvalues of Dw as well as the support of the level density of the real eigenvalues is of the 
order while the number of the additional real modes increases linearly in a. 

In the continuum limit |a| ^ 1 the scaling behavior is crucially different. Then the spectrum can be approximated 
by the microscopic spectrum of the continuum Dirac operator without zero modes plus the spectrum of a finite 
dimensional GOE. The GOE describes the broadening of the former zero modes which scale like a in the continuum 
limit and stands in contrast to the thermodynamic limit. Moreover, we have on average additional real modes. 
Hence additional real modes are strongly suppressed for larger topological indices ly as is the case for QCD with three 
colors. However the suppression is much weaker than for three color QCD. The reason is the weaker level repulsion 
which is linear for the two color case. This repulsion is also the reason why no separate peaks are visible. 

Surprisingly the continuum limit is not uniform for the microscopic level density of at the origin. In the 
continuum we would expect p^ijn = 0, A —>■ 0, a = 0) = 1. However the continuum limit yields P 5 {m = 0, A = 
0, a —)■ 0) = l/-\/2- This shows that one has to be careful when considering observables which essentially depend on 
the eigenvalue nearest to the origin for configurations with ly = 0 and also v = 1,2 if the quark mass is non-zero. 
The mass dependent chiral condensate of the quenched theory confirms this statement. It can be quite accurately 
approximated with continuum QCD without zero modes and the finite dimensional GOE for larger values of the 
index, while it has the strongest deviations for small topological indices. The ly zero modes broadened along the real 
axis push the spectrum away from the quark mass m which is the point where the limit is non-uniform. 
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As is the case for QCD with three colors, we have also found that the low energy constant W% has to be positive. 
The reason is that the chiral partition function with an even number of dynamical quarks is not positive definite if 
Wg < 0. Effects from the low energy constants Wq and WV were not considered in the present work. They may 
weaken this conclusion, though we would not expect this due to the analogy with the case of three color QCD. 

Technically the two color case is much more complicated than the three color QCD. Particularly, there are no 
analytical results for the spectral density of D^. At least for small lattice spacings the density of the real eigenvalues 
is tightly constrained by the distribution of the chirality over the real eigenvalues of Dw and by the mass dependence 
of p 5 at the origin. It is not clear if it is at all possible to derive analytical expressions for the microscopic spectral 
density of Dw. For the three color case our derivation relied heavily on the existence of an Itzykson-Zuber integral 
over a non-compact group. Such an integral is not available for two color QCD. But there is another way to deduce 
some results for the real eigenvalues of Dw in two color QCD. We expect that the results of three color QCD should 
also apply to the two color case. Mean field results generally do not depend on the Dyson index /3 d- Indeed we have 
seen that apart from factors of 2 for continuum QCD. This is sometimes referred to as orbifolding |67l. 

Let us summarize. Our results provide an analytical control of the smallest eigenvalues of the Dirac operator which 
can potentially compromise numerical simulations. It is noteworthy that the absence of level repulsion from the origin 
for the case of A^c = 2 QCD makes this effect much more pronounced than in the case of A^c = 3 QCD. 


ACKNOWLEDGMENTS 

We thank Gernot Akemann and Kim Splittorff for helpful comments. MK and SZ are financially supported by the 
Alexander-von-Humboldt Foundation. JV and SZ are supported by U.S. DOE Grant No. DE-FG-88ER40388 and MK 
is also partially financially supported by the CRC 701: Spectral Structures and Topological Methods in Mathematics 
of the Deutsche Forschungsgemeinschaft. 


Appendix A: Some Properties of Bessel Functions 


In this section we recall some properties of Bessel functions which play a crucial role in the spectral statistics of 
chiral perturbation theory of QCD. More information on Bessel functions can be found in [HH] . 

Especially, we need the Bessel function of the first kind. 


/ TT ^ 

7 ^ exp[—fxsint^Je*'^'^, for ^ G Z and a: G C, 

-TT 27r 


(Al) 


the modified Bessel function of the first kind 


/u(x) = [ ^ exp[a;cos(/j]e*'^'^, for /r G Z and a: G C, 

J-7V 27r 


and the modified Bessel function of the second kind 


Kfi{x) = 7 , / (itexp[—a;cosht]e for fj,,x G C and Rex > 0. 

^ J — OO 


The Bessel function is particularly important because of the integral 


/ OO 

dtexp[—2{e — iLx) cosht -|- 2iLy smht\e ~= 

-OO 


^ /•OO 
2 


x + y 


x-y 


m/2 




with the two angles 


(A2) 


(A3) 


(A4) 


(j)± = -L- (sign(x + y)± sign(x - y)) ■ 


(A5) 


We have assumed that £ > 0 is infinitesimally small, L = ±1 and x and y are real. The right hand side of Eq. (A4) 
can be obtained by the complex shift t —>■ t -|- ln|(x — t)/{x -I-1)|/2 — i(j)-. The angles satisfy the relations 


(j)+ = —L^sign(x)0(x^ — y'^) and e^^^~ = iLsign(x — y). 


(A6) 
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This results in the very useful identity 

[ dtejq)[—2{e — iLx) cosht+ 2iLysmht]e~^* = - -^^ (Vlv^ ~ — a;2|e®‘^+). 

2 (y-x-lLe)^^ V 7 

(A7) 

All Bessel functions satisfy a symmetry regarding the transformation fi ^ —fj,, 

J-^,{x) = = {-l)^J^ix), I-^,{x) = /^(x), and K^{x) = K_^{x). (A8) 

Except for the Bessel function this relation is only valid for integer y. 

The functions and have a particular simple representation as an absolutely convergent series which extends 
its index to general y &C, 




(-1V 


j=o 


jlT{y + j + 1) 




and 


OO 


= Y. 

j=o 


1 

jlT{y + j + 1) 




(A9) 


These series representations make it obvious that the two Bessel functions Ju and A are related by 




(AlO) 


Also the Bessel function is related to the other two Bessel functions. This relation is based on a representation of 
the Bessel functions as Meijer G-functions [68] which are given by contour integrals, 


x^J^{2x) = J 


ds r[/r - s] 2s 
27rz r[l + s] 


and 


x^K^,{2x) = \ “ s]r[-s|a;- 


^2s 


(All) 


(A12) 


The contour C encircles the positive real axis clockwise and thus all poles of the Gamma functions. When choosing 

n t. 

ds 


X = imaginary (n = ±1 and x' > 0, cf. Eq. we can take the imaginary part of Eq. (|A12|) and obtain 


Im 




= -Im 
2 

= -Im 
2 

TT 

= —n— 


r[/j, — s]r[—s]e' 


2Tri 

ds r[/j, — s] 


,27rns^/2s 


2Tri r[l + s] sin(—TTs) 
f ds r[y- s] ^, 2 ^ 


gj7rns^/2s 


2 Jq 2TTi r[l + s] 
= -n^x'^J^{2x'). 


(A13) 


We have used the reflection formula of the Gamma function r[ 2 ]r[l — z] = Tr/smirz. Since the integral is equal to a 
sum over the residues at the poles of the Gamma functions which are real, the imaginary part is obtained by replacing 
exp[i7ms] —>■ sinTrns = nsiuTrs. 

One can also understand the relation (A13) from the logarithmic cut along the imaginary axis of which follows 
from the series expansion 


1 f2Y iy-k-l)\ 




k^O 


k\ 


Z 

'T 


+ (_l)M+iin(f) j^(^) 


(-1)^ 


OO 

(i)'i: 


fc=0 


ip{k + 1) + ip{y + fc + 1) 
k\{y + ky. 


(A14) 


where ipik) = c}fclnr(fc) is the Digamma function, meaning the logarithmic derivative of the Gamma function. This 
series can indeed be calculated by taking the residues of the contour integral (A12). 

Another relation which is less known but essential to derive results for two color QGD is based on an integral over 
the Bessel functions and J^. In particular we consider the imaginary part of the integral 


/ OO 

dyy°‘KYxy) 


(A15) 
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2 

with /i, Of S K+, X G C, and Rea; > 0. By introducing an auxiliary Gaussian in the integrand we can extend 

the definition to purely imaginary x = (n = ±1 and x' > 0) and ensure that the integral over y is always 

absolutely convergent. Then we can perform the following calculation 


Im 


lo ~ / ) 

TT r[(^ + a + l)/2] 


(A16) 


Eq. jAlal l TT 

——' —n—x ^ hm 


2 e->-0 

fl 


= J dyy°'J^(2x'y)- 

The function iFi is the hypergeometric function, 


”4 r[/r+l] J™g(M+a+l)/2 


li^i 


/r + a + 1 . x'- 

---;m + 1; — 


iFi 


y, + a + 1 


-;y 


= 1 : 


+ i)/2]r[/r +1] f x‘ 


^^j!r[(/r + a + l)/2]r[/r + j + l] 


„'2 \ a 


(A17) 


The limit e —>■ 0 can be performed by rewriting the Gamma functions as an integral where the limit is trivial, 


TT r[(/r + a + l)/2] 

4 r[/i + 1] £->-0 e(/^+“+i )/2 

TT 1 x'^'^ 

4 r[(p, — a + l)/2] £->-0 e(A‘+“+i )/2 




fj, + a + 1 


; /r + 1; 



Jo 


TT r[(/r + g + l)/2] 
4r[(M-a + l)/2] 



(A18) 


The intermediate step is only true for y — a > —1. However one can analytically extend this identity to arbitrary /i 
and a since both expressions are analytical in the combination y — a. They are also bounded by the analytic function 
x'^““/r[(/x — Q!)/2] on the positive complex half-plane Re — a) > 1. Hence Carlson’s theorem, see [SI], can be 
applied which allows such a unique analytic continuation. 

Summarizing this calculation we have 


Im 



dyy^^K^ (2F-^/^x'y) 


n-x 




dyy'^Jfi. (Sx'y) 


TT r[(/x + a + l)/2] 
4r[(/r-a + l)/2] 


(A19) 


It becomes important to find the well-known results of continuum QCD from our calculations. 


Appendix B: Evaluation of the Partition Function 


The quenched partition function ( |5^ , in particular its imaginary part, plays a crucial role in the spectral statistics 
of the Wilson Dirac operator. Since its evaluation is cumbersome we split the calculation in the integration over the 
four Grassmann variables (a, a*,/3,/3*), see subsection B1 and in the explicit calculation of the imaginary part in 


subsection |B 3[ In subsection |B 2| we present some properties of the integrals involved. The explicit expressions for 
the spectral observables are summarized in subsection |B 4[ 


1. Integration over the Grassmann Variables 


We evaluate the Grassmann integrals of the supersymmetric partition function (55) using the representation of U 
given in Eq. (561. To calculate the inverse of U and its traces thereof we split U into a numerical part Uq which is a 
diagonal matrix and a nilpotent part G, 


/ 6*“^ 0 a* /3* 

I 0 -a -P 

I a a* 0 

\ /3 /3* 0 


t/o + G. 


(Bl) 
























For this purpose it is suitable to know what the square of G is, 
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{a* a + P*/3)l2 0 0 

= ( 0 2aa* a/?* + /3a* 

0 a/3* + /3a* 2/3/3* 


(B2) 


Additionally we need 


(Uo^Gf = 


iu.-^GUo^r = 


(g-*<p-si^*a + e-*‘^-'*"/3*/3)l2 0 0 

0 2e-*^-^iaa* (a/3* +/3a*) 

0 e-*‘^-'*"(a/3*+/3a*) 2e-*‘^-®"/3/3* 

(e-2*‘^-2^ia*a + e-2*‘^-2^2/?*/3)l2 0 0 

0 p-2iip-si-S2(^i:i* 


0 


2e-2*‘^-2^iaa* e-2*‘^-*i-«2 (a/3* +/3a*) |, 

s-2*‘^-*i-*2(a/3* + /3a*) 2 e- 2 *‘^- 2 ^ 2 y 3 ^* 


(Uq^G)^ = 2e-2*‘^-«i-"2 diag(l 2 , -l 2 )a*a/ 3 */ 3 , 


(Go"^GGo"^)2(Go"^G)2 = e-^^v-si-s 2 diag([e-®i + e-®=^]l 2 , -2e-®^ -2e-"^)a*a/3*/3, 
[(Go"^G)2Go"^]^ = 2e-2*v-«i-"2 diag(e-2ivi2,-e-®i-"2l2)a*a/3*/3. 


(B3) 


The inverse of U is given by a finite geometric sum 

G-i = [1 - Gp-iG + (Go-1G)2 - (Go-1G)3 + (U^^G)^] U^\ (B4) 

Now we are ready to compute the traces. 

Only even powers of the Grassmann variables and, thus, even powers of G contribute to the traces, 

Str (M + X)U = Str (M + X)Uo 

= 2{m + Xo)e*‘^ — (m' + Xi){e'^^ + e^^), 

Str(M - A)G-i = Str (M - X){Uq^ + {U^^GfUQ^ + (Go"^G)‘‘Go"^) 

= 2{m — xo)e“*‘^ — (m' — xi)(e“'*^ + e”'*^) + 2[(m — a;o)e“2*‘^“^i {jn — a;i)e“*‘^“2«i]Q,*Q, 
+2[(m — a:o)e“ 2*‘^“®2 _|_ (3 

+2[2(m - xo)e-*‘^ + (m' - xi)(e-"i + e-"=)]e-2*‘^-"i-"^a*a/3*/3, 

StrG2 = Str(Go2+ G2) (B5) 

= 262*-^ - e2®i - e2"i + 4(a*a + /3*/3), 

Str G-2 = Str + Str {U^^GU^^f + 2Str Go"2([/-1g) 2 + Str [(Go"^G)2Go"^] 2 
+2Str(Go-iGGo-i)2(Go-iG)2 + 2StrGo-2(c/-iG)‘‘ 

= 2e-2i<<2 _ e-2«i - e-2^1 

+4(e-2*v + + e-2«i)2e-i‘/2-'*ia*a + 4(e-2*v + g-*v-52 e-2«2)g-*v-s2^*^ 

+4[3 + e-2*¥’ + + 2e-*‘^-®i + 2e-^^-‘^^ + e- 2^1 + e-2*i]e-2*v-«i-^2Q,*g^^*^_ 

We also need the Sdet G+i)/2[/ which can be expanded as 


Sdet G+^'iGu = gG+i)O 0 -Oi+^ 2 )/ 2 ) [1 _ (^ + l)(e-*^-"ia*a + e-*^-"=^/3*/3) + iy{iy + l)e-2*‘^-'*i-«^a*a/3*/3] . (B 6 ) 


Collecting everything we can integrate over the Grassmann variables which only selects the highest order polynomial 
in those. To keep the notation simple we define the following abbreviations 

tzi = [(m — a;o)e“*‘^ + {m! — xi)e~'^'^], 

tz2 = iTe“*‘^“®^[(m — xo)e~^‘^ + {m' — xi)e~‘^'^], 
tzi2 = iTe“*‘^“®^“^^[ 2 (m - a:o)e“*‘^ + {m' - xi)(e“®^ + e“®^)], 
d, = (12 + l)e-^‘P-^\ 
d2 = (12 + l)e-*‘^-*^ 
di2 = K^^ + l)e-2*‘^-*l-^^ 

tai = 4a2[l + g-*¥>-^i(g-2*v + g- 2.1 g-*¥>-«i)], 

t„2 = 4a2[l + g-*¥>-^2(g-2^v, g- 2.2 g-i¥>-«2)]^ 

t„i2 = 4a2e-*2‘P-si-s2(3g-2i¥> g-si-s2 g- 2.1 ^- 2*2 + 26"*'^-"^). (B7) 
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The final result for the generating function is up to an overall constant 


Z oc 


s/:-. 


dsi 


sinh 


Si - S2 


„(j/+ 2 ) 2 ¥> 2 )(si+S 2)/2 gjj-p _j_ 


X exp \—2iLxQ cos (p + iLxi{coshsi + coshs 2 ) + 4a^ cos2i^ — 2a^(cosh2si + cosh2s2)] 

X [(c?i — tzl — tal)id 2 — tz2 — ta2) + ^12 — <^ 1^2 + tzl2 + fal2]- 


(B 8 ) 


This expression is explicitly shown in Eq. (631 which immediately makes clear that the integral splits into a finite sum 
and each term factorizes into a compact one-fold integral and a coupled, non-compact two-fold integral. 


2. Properties of the Functions 5^,a and <1>^ 


In this appendix we discuss properties of the integrals and 4)^ which were introduced in Eqs. (651 and (64|. 
Their properties are in particular useful for the numerical evaluation of the integrals. 

The integrals (64) and ([ 6 ^ satisfy the following recurrence relations 


0 = -f (xo -h to)4>^+i -f {xq - - 4a^$^+2 -h 4a^$^_2, 

0 = + (xi + m')S^-i^a+l + (a^l — + i — 4a^(2iS'^_2,a-|-2 ~ ‘^S^^2,a+2 — S^-2,a + S^+2,a)- (B9) 


These relations can be found by acting d^p and {dg^ +^ 82 ) onto the integrands which vanish under the integrals. They 
are useful to check the quality of the numerical evaluation of the integrals (64) and (65). 

Two other useful relations of S'^,a concern the derivatives which are important to calculate the observables in¬ 


troduced in subsection 

by 


IIB The derivatives with respect to the quark mass m' and the axial mass xi are given 


Moreover, there are also relations between a negative and positive index /i, 

^p{m,Xo,a) = (-l)'"$_^(-m,Xo,a) and S'^,a(w', Xi, a) = (-I)'"S'_^_a(-TO', Xi, a). 
These relations are helpful when computing these integrals. 


(BIO) 


(Bll) 


3. Real and Imaginary Part of Sp^o 


We assume /i > 0 because of the relation (Bll). This simplifies the computation of the imaginary part of the 
integral (65). This integral can first of all be written in a more suitable version. For this purpose we introduce a 
Gaussian integral over an auxiliary variable t to linearize the term cosh 2s = 2sinh^ s -I- 1 to a sinhs term. 


Sp^airn\xi,a) = 


i^Lr r 

VStto^ J-c 
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dte->^^ 


,, q ;+1 


^ 22/2 _ 1 


X exp 


8a2(22/2 _ 1) 


(f — m')’^ + 2iLty sinh s -I- 2iLxiy cosh s — 4a^(2y^ — 1) 


(B12) 


Hence the variable t acts as an effective quark mass. The integral over s is equal to a modified Bessel function of the 
second kind. 


Sp,^aim',xi,a) = 
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dy 


dt 
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a+l 


1 


TTO^ 


72?/2 - 1 {t-xiY 




X exp 


8a2(2?/2 - 1) 


(t - m')^ - 47(27 - 1) 


K, [2y^\t^-x\W<^'^ (B13) 


with (j) = —Lsign(xi) 0 (x 2 — t^) 7 r/ 2 . The phase is important and reflects the transformation (only a complex shift in 
the variable s) to bring the integral into the form (A3). 
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The imaginary part of this integral consists of two contributions. The imaginary part has to be taken either of the 
term l/(t — xi)^ which results in the (/i — l)st derivative of the Dirac delta function, — xi), or of the term 

which can be dealt by relation (A13). Thus, we find, 


lim Im S'u afm', Xi + ie, a) = (—1)'^ ^ 

e ->-0 



dy 


/•OO 2/“^^ 

/ dt — , exp 

^oo 


8a2(2?/2 - 1) 


(t — m')^ — 4a^(2j/^ — 1) 


J_ ^ (/r — fc — 1) ! 


[y‘^{xl - - a;i) - signxi (^2y^ xl - 


Xi+t 
xi — t 


li/2 


Q{x1 — t^) 


(B14) 


The sum in the first term comes from the series representation of see Eq. (A14|, since the {fj, — l)st derivative of 
the remaining parts vanish at t = xi. 

The real part of the partition function Z is needed to compute the chiral condensate E(to). The real part of the 
expression ( |B13 ) cannot be easily taken since the integrand has a pole of order y at t — 0. Therefore we integrate by 
parts and symmetrize the integrand with respect to t —>■ —t. Then the real part is given by 


Re Xi = 0,a) = ,xi = 0,a) 

1 


dy 


dt- 


\/2'Ka?{y — 1)! Jo y/2j/2 — 1 


:exp[-4a2(2y^-l)] 


(B15) 


X -djr 

t * 


\t' + 1\^ exp 


'8a^(2y^ - 1) 


(t' + t — m)"^ 


{2y\t' + t\) — {t —> —t} 


J t'=o 


Note that </) = 0 at xi = 0. Despite the modulus in \t' + (2y\t' +1\) the function is differentiable at t' +1 = 0 

because the pole of the Bessel function is cancelled and the derivative of the logarithm yields a pole which cancels 
with the zero of the Bessel function cf. Eq. (A14). 


4. Explicit Expressions for Spectral Observables 


In this appendix we give explicit expressions for the spectral observables in terms of the functions and 


which were introduced in Eqs. (65) and (64), and we use the short-hand notations 




= S^^aijn, A = 0, a), and A, a). 


(B16) 


For the spectral density of the Hermitian Wilson Dirac operator we hnd 


P 5 (m, A, a) = 5a;Jm2'^(M, A,a) 


M=m]L4 


X^(A+ie)l4-J-Al4 

= 160“^ [$i,_4(5,y+i4 — iSy_iq)-l-2$i,_3(5y+2,2 — 5^,2)+‘hi/-2(45^+3,3 — ASv+1,0 — ‘5iy+3,i -I- iSiv+iq) 

+ (iS^+4^2 — Su+ 2,2 + Su ,2 — S,j-2,2)+^i^{Siy+5,l — ‘5i/+ 3P -|- 4iSy+i_3 ~ 4iSi,_i^3 ~ 2S^+iy -|- 2iS^_i,i) 

+2$j^_i_i (4iSy_i_2_4 — 4:S^^4 — 8S^+2,2 + 3iSy^2)+4’;y+2(‘5i/_i4 — 5^_3p)] 

+4a^ [$j^_2((2j 2 -|- l)5y_|_i^i — (2i2 — l)5y_i4)-l-2$j^_i(5y_|_2_o + (12 — l)5i,+2,2 — {i' — l)5i,_2) 

+4’[/((2z 2 -|- l)5,y+in ~ {2v — l)5i,+3q + 4i2iSj/+3^3 — 4z/iS;/+i^3)-l-2$^4_i((i2 -|- 3)iSi/^2 — {v + l)5i/_2,2 ~ 

—8aJ{m — A) [<i)y_3 — 5iy_iq)-l-2<i)^_2 (5y+2,2 — ‘5iy^2)+2$j^_i(5i/+3_3 — 5i/+i_3) 

-\-^v(Sy^2 — '5j/_ 2,2 + '5i/+4^2 ~ '5y-|-2,2)+4’jy+l (2iSy_|_i_3 ~ 2iSi/_i_3 ~ -f 5^_iq)] 

+ (m — A)^ [$y_2(5i,+iq — S^-iy^+ 2 ^v-i{S^+ 2,2 — ‘5i/,2)+4>jy(5jy+3q — 

—2(m — A) -|- l)iSi/+iq — +^,^(^3,^+2 ,0 + vSi, 4 . 2,2 — i^Si,^ 2 )\ 

+ 3)5i/+iq — {v + (Blf) 
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The analytical result for the chiral condensate for the quenched theory is given by 


I](m, a) 


dm'Re Zi,{M, X, a) 




— 16a^ +2<i>i/_3 (s' i/+2,2 + S,y^2 )+$^_2(4S + 4Si/+i^3 ~ S 1/ + 34 ~ Sjy+lj) 

+24>i,_i(S^+ 4^2 + Si/+2,2 + Sjy^2 + Si^-2,2) + 4>^(Sjy4-54 + Si/+3^i + 4S^+i^3 + 4Si/_i_3 ~ 2Sjy+i,i — 2Sjy_i_i) 
+24 >i,4_i( 4S^4-2,4 + 4Sy4 — 3Sjy+2,2 ~ 3Sjy^2) + 4>,y4-2(Si/_i4 + S;/_34)] 

+4a^ [<i>y_2 ((2 j/ + l)S,y+i4 + ~ l)Sy_i^i)+2$i,_i(Si/4-2,o + {v — l)Sj/+2,2 + {i' — l)Sy_2) 

+ $^( 4 l^Sy+ 3_3 + 4 i/S;^ 4 _i ^3 — ( 2 l/ — l)Sy4_3j — ( 2 ^ + 1 )S^ 4 _i,i)+ 2 $j^_|_i (^{1/ + 3 )Si ^^2 + + 1 )S^_ 2^2 ~ ' 5 'i/,o)] 

—8a^m [$y_3 (Sy_i_i 1 + S,y_i^i)+2$[^_2 (S i/+2,2 + )+24)y_l (Sj^4_3^3 + Sjy4_i_3) 

+ 4*y(Sy^2 + Si,_2,2 + S;/4-4^2 + Sj/+2,2)+4>y+l(2Sy+1^3 + 2Sy_i^3 ~ >S'j/+l,l ~ Sy_i^i)] 

+m^ [$i._ 2 (S^ +14 + Sjy-i^i )+2$y_i (Si/+ 2,2 + S^_2)+4?jy(Si/+3^i + S 1 /+ 14 )] 

+ l)Sy4-i4 + z^Si/_i4)+$jy(Si/+2,o + i^Sjy+2,2 + J^Sjy_2)] 

+v^^{{v + 3)S;^4_i^i + + l)Si,_i.i). (B18) 

The distribution of chirality over the real eigenvalues can be written as 


Px(w,a) 


-dmha Zy{M, X ,a) 

TT 


M=ml4 


X^{A+ie)l4-^Al4 

= Ida'^ H " 2 $ i /_3 ( tSi ^+ 2,2 + ^ 1^,2) 1^-2 ( 4 tS iy+3,3 + 4 S ^+ i _3 ~ ' 5 i /+34 ~ ‘ 5 ^+ 1 ,!) 

+2‘l>y_i(Si^4_4^2 + ‘5;y+2,2 + <Su,2 + S^- 2 , 2 ) + 4>jy(Sy4_5_i + S^+s^i + 4S^4_i^3 + 4Si^_i_3 ~ 2Si^4_i^i — 2Sy_i,i) 

+ 2$y+i (45^4-2,4 + 4 Si/4 ~ 3Si/ 4-2,2 ~ 3S,y^2) + 4>y4-2(Sy_i,i + Si,_34)] 

+ 4 a^ [$y_2 (( 2 i/ + l)Siy+i,i + { 2 v — l)Si/_i4)+2$j^_i(Si/+2,o + {v — l)S^+2,2 + {i' — 1)S^^2) 

+$j^ (4i^Sy+3_3 + 4:^Sy4-i_3 — (2 z 4 — l)Si/+3^i — {2iy + l)Sy+i4)+2$j^_|_i ((z/ + 3)Siy^2 + (y + l)Siy_2,2 ~ 

— 8a^m[$i^_3(Sy4-l^l + S^_14)+2 $i^_2(Si,+ 2,2 + S^_2)+24>y_i (Sy4-3^3 +Sy4-1^3) 

+4^!^ (Sy ,2 + Sy_ 2,2 + S^+4^2 + ‘5i/+ 2,2 )+4>y+4 ( 2 Si/_|_i ^3 T 2 Sjy_i ^3 ~ ~ 1 , 1 )] 

+771^ [<i)j ^_2 (Sy +14 + Sy_l^l)+ 2 $J ,_1 (Sy + 2,2 + Sy^2 )+4’i/(S;^4_34 +Si^4_1^i)] 

- 2 m[<^v-i[(y + l)Si. + 14 + l^Sy_i,i)+$y (Sjy+2,0 + i^Sy+2,2 + vSv^ 2 Y\ 

+z^$^((z^ + 3 )Sj/4-i^i + [v + l)Sj/_i^i). (B19) 

Though all three expression look quite complicated they are a finite sum of two kinds of integrals, only. This simplifies 
the numerical evaluation a lot. 


Appendix C: The Continuum Limit a —> 0 


To be self-consistent we briefly review the exact continuum limit (see subsection C 1) and the Gaussian orthogonal 
random matrix ensemble of finite matrix size (see subsection C 21 which should describe the spectral broadening of 
the former zero modes into the real axis quite well at small lattice spacing |a| <§; 1. In particular we wish to show 
how to extract the known results for the spectral observables from our calculations. 


1. Exact limit a —0 


In the continuum limit a —> 0 the partition function takes the form 

+ (M, X, a = 0) = (to - xo)^4)y_2+,o + 2 (to - xo){m' - Xi)^^-iS^+i^i - 2v(m - a;o)4)y_i+,o 
+(m' - xi)^$y+ 4 _ 2 ,o - 2u{m' - -I- (z^ -b l)z7$;,+_o- 

with 

m/2 


$;_,(TO,xo,a) = 


m — xq 
m xq 




(Cl) 

(C2) 
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and 


= 2 —-— / 


dyy° 


'’'2 — K^[ 2 yJ\'m''^ — x\\e 


t/2 — ^21 


(C 3 ) 


with (f) = —Lsign(xi) 0 (xf — m'^)7r/2. We have to be careful about the phase e*'^ since it may change the sign of the 
result while it is unimportant for the phase of yjrn? — Xq in the compact, analytic integral. The reason is the cut 
along the negative real line of The index a only takes the values a = 0,1 which simplifies the calculation a lot. 
In the first step we simplify the partition function by employing the relations 


+ (a;o + + (a;o - = 0, 

zK^+i{z) + {zd^ - y)K^{z) = 0 , 

K^+^{z) + K^_i(z) + 2d,K^{z) = 0. 

The second and third relation is needed for expressing S';y+i,i and in terms of 8 ,^^ and Ki, respectively, 

(v4 


(C4) 


(to' - xi)S^+i^i{m',xi,a = 0 ) = ^ ^ ^ ,xi,a = 0 ) + 

(to' — xi)^5'^+2,o(m.', a;i, a = 0) = —(to'^ — x\)S^fi{m', xi,a = 0) 

, v+l 


|to' 2 — xf |e*‘' 
(to' — Xi)" 


-K, ( 2 A/|TO' 2 -x 2 |e'^ 


+2 


to' 2 — x?|e"' 


K,+A2J\m'^-xlW^ 


(to' — Xi)" 

resulting from integration by parts. Then the quenched partition function reads 

^-\/|to' 2 — xf (to — xo)"/^ 


(C 5 ) 


ZAM,X,a = 0) = 2- 


(to' — xiY{m + xo)''/2 


[(to - Xq) - (to - xj)]/^ 2 Wto2 - 


X / dyK^ 2 yJ |to '2 - xf |e"P + \ m? - x§ 4 +i 2 Jw? - x§ [ 2 J |to '2 - xf je* 


+ -y/|TO'2 - xf|e"^/y (^^2^Jw? - Xg'j K^+i ^2y^|TO'2 - xfje"^ 


(C6) 

The last two terms also appear in QCD with three colors. In this way one can easily check the normalization 
Z,^(to 14, XI4, a = 0 ) = I. 

In the second step we take the derivative with respect to to' or xi and set to = to' and xq = Xi = x, yielding 


dm'ZAM,X,a = Q) 


M—7n\i 

X=x1a, 


— Amli. 


2m 


dy [2V|TO2-x2|e'^4+i ( 2/1 


■\/|to 2 — x2|e' 

- 4 to [/,+! ( 2 V|m 2 -x 2 |e"^) ( 2 V|m 2 - x 2 |e'‘^) + A (2/I 


— 4 m/, 


X — m 


.(2v^ 


^2-\/|to 2 — x2|e''^^ J dyKi, ^2j/-\/|to 2 — x2|e''^^ 

14 (j/) + ylv (2a/|to 2 - x2|e'‘^) (y) 

f dyK„ {^yVl 


m? — K, 


m? — x‘^\e*'^ K, 


( 2 ^ 


t/—2-y^|m2— 


(C 7 ) 


and 

d^,ZAM,X,a = 0) 

2 x 


vm 


M—7n\i 

X=xl4. 


■ 4 x 1 , 


(2V|m2-x2|e''^) 


dyK^ 


0 


+ 4 x 


■\/|to 2 — x2|e''f' 

{ 2 -/\ 


dy 


2v'|to 2 - x2|e"^/i.+i (^2 v'|to 2 - x2|e'‘^^ K,, (y) + yl^ ( 2 y/\m‘^ - x2|e'‘^^ (y) 

(2y[ 


4 x/„ 


11/+1 


TO — X 

to 2 — x2|e"^ 1 itTiz-i ^2 '\/|to 2 — x2|e'‘^ ) + /, 


( 2 ^ 


dyKi, 


m? — x2|e*'^ 1 K, 


( 2 ^ 


J y— 2 ^J\ra^ 


(C8) 
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Up to an overall sign, Equation (C7) is equal to the chiral condensate E when setting x = 0 (implying ^ = 0), i.e. 
the chiral condensate without zero modes is 


^chGOEM =4m 


I^+i{2\m\)K^_i{2\m\) + I^{2\m\)K^{2\m\)+I„{2\m\) / dyKi,{2y\m\) 


(C9) 


This result is in agreement with the expression in |66j , but has the advantage that there is no need to distinguish even 
and odd v. For odd v the Bessel function K^, can be written as a total derivative allowing us to evaluate the integral 
exactly, 

K2k+i{x) = -2^[K2kix) - K2k-2{x) + • • • + i-lfK^ix) + {-l)'‘+^lKoix)]. (CIO) 

ax 2 

For u = 2k this results in the condensate 

E(m) = -- + AI,{2m)[K2k{2m) - X2fe-2(2m) + • • • + {-l fK2{2m)] + (-l)'=+n,(2m)Ko(2m) 
m 

—4mIy+i{2ra)K^-i{2m) — AmIi^+i{2m)Ku-i{2m), (Cll) 


which up to a rescaling m —>■ 2m and E —>■ E/2 agrees with the result in 
relation for Bessel functions but we are left with an integral over Kq. 


1. For even ly we can also use the recursion 


K2k{x) = 2{-ir/^ — 


'v/ 2-1 


fc=0 




(C12) 


This results in 


E(77i) 


- 477i/i.+i(2m)iFiy_i(2m) — Amy+i{2m)K,y-i{2m) 

m 


U2-1 -oo 

+4(-l)''/^/^(2m) F] (-l)''-?^ 2 fe+i(a;) - 2TO(-l)''/^/,.(2m) / dyKo{2my). 
k=o 


(C13) 


The integral over Kq can be expressed into modified Struve functions, and it can be numerically shown that it agrees 
with Eq. (8) of [5S] . 

Let us underline that the valence quark mass dependence of the chiral condensate for chCSE is also de¬ 
rived in |66j . The authors again obtain separate expressions for even and odd v which can be simplified to 
- 2 xK]srj+ 2 ui‘^x) dylN{+ 2 ^i‘ 2 xy) up to rescaling. 

From the result ( |C7[ ) we can immediately deduce the well-known result for the distribution of the chiralities over 
the real eigenvalues which is a Dirac delta function in the continuum limit. 


Pxi'm,a = 0) = dm'a = 0) =v5(m). 

X =—2e]L4—>-0 


(C14) 


Only the first term of Eq. (C7) contributes to this results since the remaining parts are real for m > x = 0. 

The spectral density of the Hermitian Dirac operator Dq can be obtained by setting x = \ + le and taking the 
imaginary part of Eq. (C8) in the limit e —)■ 0. We have to distinguish two cases, |m| > | x| or | m| < |x|. I n the 
first case all terms are real and vanish while in the latter case we have to apply the relations (A8|, (A13l, and (A19) 
yielding, 


Pq{\, m, a = 0) = —Im X,a = 0) 

TT 


AI—7nl.4 

X^{X + i€)l4^Xl4 


= i'S{m — A) + 


|A| 

yx^ — m'- 


:JX2VX^-m^) 


p2\/ — m? 


1 - 


dyJv{y) 


0(|A|-|m|) 


-l-2|A| j/(2\/A^ — m^) — Jy+i(2-v/A^ — m^)Ju-i{2\/ A^ — m^) 0(|A| — |to|). 


(C15) 


cf. [SS]. Note that Eqs. (A13) and (A19) are multiplied with n 
TO'^)7r/2 for L = 1. 


—sign(A) because of the phase ^ = —sign(xi)0(xf — 
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2. Spectral Observables of the Gaussian Orthogonal Ensemble 


We choose the abbreviation Y = (M — X)/{Aa) = diag(2/o, 2/o, 2/i) 2/i)- The partition function (62) is equal to a 
partition function of x real symmetric random matrices distributed by a Gaussian. Let Sym(^) be the set of these 
real symmetric matrices. In terms of these matrices, the partition function reads, 


4';^e (t) 


1 


i)/47ri^/2 

For the cases ly = 0,1 they take particular simple forms, 


7' 


(0) 

GOE 


(?) 


Met(iL-yil,) 


= 1 


and 




( 1 ) 

GOE 


(?) 


= 4 / dE^——e ^ = 1 + V^ivo - - j/o)e 

Vtt E-yi i 


(C16) 


(C17) 


(CIS) 


In the latter equation we assumed that the imaginary part of xi = Rexi +iLe {L = ±1) is infinitesimal small, e —>■ 0. 
The error function erf(z) = 2 dz'e~^ jspK is real despite the imaginary unit. 

For u > 1 the result directly following from the representation (C16l are not that trivial as the ones in Eqs. (CI7) 
and (CIS). Therefore we start for those cases from the supersymmeric representation (62) which, in the parametriza- 
tion (59), reads as 


7(<^) 

^GOE 


( i nOO pOO pOO 

Y]=— du dvi / du 2 |ui — +'"1+^2) — 
7 tt (j/i 


( 2/0 - iuY 


1-^ 


2?/i -V 1 -V 2 


+ 


- VlY^'^iVl - V2Y^'^ 

v{v-l) I 


4 (yo - *u)(yi - vi)(yi - ^ 2 ) 16 {yo - iuYijJi - vi){yi - V 2 )_ ' 


(CI9) 


We have already integrated over the orthogonal matrix O and the four Grassmann variables 77 , y* 1 X^ X* ■ Moreover, we 
rescaled the remaining variables (u, Vi,V 2 ) —>■ 4a(M,Ui,U 2 ). The normalization constant is fixed by yg = yi = ie ^ 00 . 
The integral over u is equal to Hermite polynomials, 


due 


-2u‘ 


(yo - 444)^ = J-2-^77^ (2yo), 


where we have chosen the monic normalization Hf^ (z) = + ... Then the partition function is equal to 


^GOE (^) = diu (2yo)I,.,o(yi) + 77,x-i (2yo)Ii.,i(yi) + iy{iy - l)i7j.-2 (2yo)X,,+ 2 .o(yi) 

dv2\vi — . 


with 


T^.a(yi) = ^7=9° 




dvi 


(yi - 'ci)^/^(yi - n 2 )^ 72 - 


(C20) 


(C21) 


(C22) 


The latter twofold integral exits due to the non-vanishing imaginary part of yi = m — Rexi — iLe. This integral can 
be brought into a form which is easily integrable, and the spectral observables we are interested in can be readily 
extracted. To find such a representation it is suitable to understand the integral (C22) in terms of an integral over a 
2x2 real symmetric matrix H, 


T t \ aa f ^dle 

7sym(2) det'^/ (yi 


u-HY 


(C23) 


where vi and V 2 are the eigenvalues of H. The determinant in Eq. ( |C23 ) can be written as an integral over 2x2 
positive definite real symmetric matrix H' 


dee/Yyi^2-H) 


JLTTfl/2 


/ exp[—2zLtri7^(yil2 — i7)]det^^ 

7Sym^(2) 


(C24) 
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The index /i has to be larger than 1 and thus v >2 which is okay since we know explicit results for ^ = 0,1. The phase 
in front of this integral is important as well. The integral over H yields a Gaussian integral in H'. After diagonalizing 
H' = O'diag{Ei, £ 2 ) 0 '^ with O' G 0(2) and Ei,E 2 > 0 we find 


^n,a{yi) — 


a)/2 roo noo 

- -—- / dAi / dE2\E^-E2\{E^E2)^>^-^^/^{E^ + E2reM-{El+El)-2iLy^{Ei+E2)\. 

[p, — Z)\ Jq Jq 

(C25) 

Already at this point we can set the imaginary part of yi exactly to zero. We order the eigenvalues Ei > E 2 and 
change to polar coordinates Ei = rcosp and E 2 = rsin(^ with r > 0 and 0 < p < 7r/4, 


^n,a{yi) — 


2(M+3a)/2gZL7r(/j-a)/2 


X exp 


(m- 2)! 
-i?2 - 


J dR J sin^^“^^/^(2i^) cos yip + — j cos“ — — j 


VSiLyiRcos ~ 

di^sin^^“^^'^^(2(^) cos (^p + cos“ 

r ('i^, 1; -2,; (^ - ^)) 


2(/J.+3a-2)/2 ptt/A 

(m-2)! i 


+ V8r - 0 M 

We have employed Kummer’s confluent hypergeometric function |68) 


^ Q; -|- 1 3 


-\-2ylcos^ [p- 0 


(C26) 


M {a,b;z) = ^ 


j=o 


r[Q + j]r[b] 

r[a]r[& + j] j! 


(C27) 


We underline that only one of the two terms in Eq. (C261 is real and the other one imaginary depending on the parity 
of p. In random matrix theory this subtle difference between even and odd matrix dimension for real matrices is well 
known |61] . 

Let us define the abbreviations 


dpsm^'^ ^^/^(2p) cos + 0 cos“ 0 - 0 M 


y + a 1 „ 2 2 

^—> 2 ;- 22/1 cos 




(C28) 


and 


^Xl{yi)=yiJ^ dpsm^>^ ^)/2(2(^)cos 0+0 cos“+0(^-0 M ^^'^i-0^0-2j/2cos0(^-0^ . (C29) 

Then the real part of the first derivative in yi and setting yi —>■ m/da yields the contribution of the chiral condensate 
which results from this finite dimensional Gaussian orthogonal ensemble, i.e. 


n(i/>i) (m 


■^GOE 


I = - dyn'ReZ, 


O) 

GOE 


(-) 


Y =:(m/4a) I 4 


(C30) 


- (t) K (£) -s (£) - (£) ^ (£)-(- m (£) 


(zz-2) 

2(‘'-i)/2r[(i/ + 3)/2] . /TTI/ 

•-;--- sm — 

{iz-2)la V 2 


12+1 


H, 




/rn\ ( 2 ) (m 

\2a) Ua 


-H, 


v-2 


/TO\ ( 1 ) zm\ 

\2a) 


Note that only one of the two terms contribute because either cos (pKv 12) or sin {'kv/2) does not vanish depending on 
the parity of v. The “chiral condensate” for the case 1 / = 0 vanishes, 




while the one for = 1 can be read off from Eq. (|C18|), 


,(,2=1) /TO 


^GOE 


V4a 


TT erf(im/4a) 


Iba^ 


' exp - 


TO 

16a2 


■ 


(G31) 


(G32) 
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Other representations of the partition function and the chiral condensate in terms of Hermite polynomials and its 
Cauchy transform exist in the literature [61]. However we chose the representation Eq. (C30l which can be easily 
evaluated numerically and the cases even and odd v can be discussed on the same footing. 

The level density is obtained from the imaginary part of the first derivative in yi, 


0>i) (^\ 


1 


Pgoe 


-j = a,„.imz, 


GOE 


2''/‘^T\(v + 2)/2] . fT:v\ 
■ sin ' 


{y — 2)!a 


V 2 y L 




(y) „ 

V } V—(m/4a)]l. 

(S) "S (£)+(S) "S (S) - (i) "3- (£) 


(C33) 


2('^-i)/2r[(^ + 3)/2l (-Kv 
- cos ' 


[y — 2)!a 


V 2 


y + 1 


-H„ ^ 


(S) "S oe) "S o 


The exceptional cases for ^ = 0,1 are 


^goe ( 4 ^) - 0 - 


and 


y=i) 

Pgoe 


(-) 

\AaJ 


VI 6 


: exp 


iTTO^ 


16a^ J 


(C34) 


(C35) 


These expressions are a good approximation for the level density of the real eigenvalues of the non-Hermitian Wilson 
Dirac operator at very small lattice spacing |a| <?; 1 which are the former zero modes of the continuum QCD Dirac 
operator. 


Appendix D: Computation of the Saddle Point Solutions (821 and (841 

We have to solve the rational equation 

p{v) = jiy- + ^(.y + y~^) - ‘^a^{y‘^ - y~^) = o. 


(Dl) 


with real coefficients. In particular each eigenvalue z oi V satisfies exactly the same equation. Because the rational 
function p{z) has real coefficients it has either no, two or four real zero points satisfying Eq. (Dl). All other zero 
points are complex conjugate pairs. Since p{z —>■ 0) —>■ 2a^z“^ —>■ +00 while p{z —>■ ± 00 ) —>■ — 2(?? —?> —00 we have 

at least two real zero points. The question is: What happens with the other two zero points? We are in particular 
interested the complex solutions since they are important for the level density when we have the general situation 
with A 7 ^ 0. 

Splitting the eigenvalue z oi V into a radial and angular part, i.e. z = re*''’ with r > 0 and tp s] — tt, 7 r[, we obtain 
two equations, 


m + A 


-r cos if 


m — \ 


r cos Lp — 2a^ {r^ — r )(2 cos*^ — 1) = 0, 


m + A m — X 


r — 4a"^(r*^+r ) cos(/J ) sin:/? = 0, 


(D2) 


corresponding to the real and imaginary part of the saddle point equation, respectively. The second equality is either 
solved by sim^ = 0 implying a real solution or the cosine can be expressed in terms of the radius r > 0 . 

Let us assume sim^ 7 ^ 0. Then the equation solved by the radius is 


sinh^ 2d+ [1- 


(8a2)2 


. , 2mA 

smh 2-i9 + -—= 0 


(8a2)2 


with r = . This equation has always one real solution which is 


sinh2iy = — 


mX 

(8a^ 


r2A2 

WvJ 


1 


(8a2)4 27 V (8a2) 


r2 — A2 


1/3 


- 1 


mX 

(8a^ 


T-2A2 

■>2\4 


1 


(8a2)4 27 V (8a2) 


z2-A2 


- 1 


(D3) 


1/3 


(D4) 
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where we have to take the following root = sign(x)|x|^/^ if x is real. There are two additional real solutions in 
the region with 


(8a2)4 ^ 27 { (8a2)2 ^ 8a‘^ ^ 


/ TO n2/3 


1 3 


(D5) 


which implies < m? — (8o^)^ and can only appear if the quark mass satisfies |to| > 8a^. However in this region the 
cosine of the angle cp has to be larger than 1 implying that it has to be complex which should not be the case. This 
can be seen in the following short computation. 

The other two solutions are given by 


sinh2d„ = —sign(TOA) 


_j_g—■i2n7r/3 


1/3 


e-»»/3 I ^ + ,sig„(„A),/ 1 (_ iV _ 


(8a2)2 


27 V (8a2)2 


toAI , , /l /to2 — A2 A"* to2A2 

-zsign(TOA)W — [ ...... - 1 I - 


(8a2H 
1/31 


(8a2) 


27 V (8a2)2 


(8a2)4 


(D6) 


with n = 1,2. Hence each solution is bounded as 


1 im^ — . , A , A . 1 „ n / to2 — A2 1 /to2 — A2 jm^ — 

— 1 < sign(TOA) sinh2di < ^-1 ^ ^—h 2 < cosh2^1 < 


ysV 


(8a2)2 


73 V 


(8a2)2 ’ 


1 /^2 — \2 Xjpi _ y2 

0 < sign(TOA) sinh 2^2 < —i= -i / ,^9 -1 1 < cosh 2i?o < —f; -i / tao —h 2. 


73 V 

We combine the two equalities in Eq. ( |D2[ ) such that we have for the angle 

(to + A)^r2 — (to — A)^r“2 


I cos2(/3| = 


73 V 

2|toA| 


32a4(r4 — r~^) 

This allows us to make the following estimations 

w? + A^ 2|toA| 


n = 1 :| cos 2ip\ > 
n = 2 :| cos 2ip\ > 


(8a2)2 cosh2i?„ (8a2)2| sinh 2d„| 

^H + |a|7 


(D7) 


(D8) 


> 


— X^ 8a?y^w? — X^ — (8a2)2 — X^ ’ 

V3(rn^ + A^) 712|toA| 73(|to| + |A|)=^ 


> 


— A2 + 2(8a2)2 — X^ — (8a^)^ 8a‘^y^m^ — A2 + 2(8a2)2 


(D9) 


In both cases the right hand side is a function which is monotonously increasing in |A| implying that |A| = 0 is the 
minimum. The remaining function in |to| is also monotonously increasing such that we have a minimum at the values 
at ItoI = 8a^ which results in 


I cos2/9| > 1 


(DIO) 


meaning that those solutions for sinh 2d are forbidden. 
Hence we have only a complex solution if 


M 

8a2 


> 


8ay 


-1 3/2 


which is 


Uq = iL with ip = arccos 


TO cosh d + A sinh d 


8a? cosh 2d 


(Dll) 


(D12) 


and d is given by Eq. (|D4j). The sign L is fixed because of the non-compact double integral over Si and S 2 - The 
infinitesimal imaginary shift iLe in A and the singularity at , e®^ =0 prevents to shift the contour to the solution 
z = in the thermodynamic limit. 
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The region (Dll) implies a spectral gap of in the interval 

3/2 


A e 


(|m|2/3 - (8a^)2/3) ,(|m|2/3-(8a2)2/3) 


3/2' 


(D13) 


cf. Fig. I^a). It is closed if \m\ < 8 a^ such that we enter a new phase which is the Aoki phase for lattice QCD with 
two colors. 

The situation of the spectrum looks slightly different at A = 0. In this case the saddle point equation takes the 
simple form p{z) = 2a^(z — z~^){z + z~^ — m/(4a^)) = 0 where we can readily read off the solutions. Plugging these 
solutions into the exponent (801 we have to choose the maximum which is 


Un 


8a2 


= i + \/l - (£ 2 ) ’ l™l < 8 a" 


(D14) 


zLsignTOl4, 


Iml > 8 a^. 


Again we would obtain always two solutions but only one is accessible in the thermodynamic limit due to the 
infinitesimal increment iLe. 


Appendix E: Density of p^(m — 1,A = 0,a —>■ 0) 


In this appendix we consider the continuum limit of p^{m = 0, A = 0,a) and show that the convergence to the 
continuum limit is non-uniform. The level density reduces in this situation to 

p^{m = 0, A = 0, a) = 16a^ [2<i>4iSip — 2$3 { 82,2 — 50 , 2 )+44>2( 53^3 — 5i^3) 

-|-2<i>i (452,4 ~ 45o,4 ~ 52,2 + 25o,2 ~ 54 , 2 )+4>o( 55,1 — 53,1 -I- 85i,3 — 45i,i)] 

-l-4a^ [2$i (25 o, 2 — 5o,o ~ 52,o)+4’o(5i,i -I- 53 , 1 )]. (El) 

Thereby we used the relations ( Bll[ ) between positive and negative indices. Considering Eq. (64), the compact 
integrals vanish for p, an odd integer if the quark mass and the axial mass are zero. Hence we have 

p^{m = 0, A = 0, a) = 16a^ [2$45i,i -I- 4<i>2(53,3 — 5i,3) -I- 4>o(55,i — 53,1 -I- 85i,3 — 45i,i)] -|-4a^<i>o(5i,i -I- 53 , 1 ). (E2) 

Thus we have only non-compact integrals 5^,^ with the index a = 1, 3 and p odd and positive. The integral over y, 
see Eq. (65), can be integrated for both indices a = 1,3, 




and 


' 5^,3 ~ Su.A — 


41,1 


27r(8a2) 
(_ 1 )( m - 1)/2 


ds 


o-MS 


cosh 2 s 


exp[—4a^ cosh 2s], 


ds- 




27r(8a2)2 cosh^ 2s 


exp[—4a^ cosh 2s]. 


We can omit the exponential in the limit a —>■ 0 if /i is small enough, i.e. 


5, 


4 t,l 


and 


8fj.,3 — 5, 


|q|<i (—1)(a‘-i)/2 j 

27r(8a2) / 

|a|<l (— 1 )(m- 1)/2 
27r(8a2)2 


ds 




4 t,l 


cosh2s 32a2 cos(/r7r/4) 

■MS (_X)(m-1)/2^ 


for 


< 2 , 


ds 


cosh^2s 8(8a2)2sin(^7r/4) 


, for \p\ < 4. 


(E3) 


(E4) 


(E5) 


(E6) 


For larger p the exponential becomes crucial since it guarantees the integrability. Before performing the other integrals 
let us see what remains to be calculated, 


/95(to = 0, A = 0, a) 


|a|<l 


16a 55,1 — 53,1 


1 


2ii/2a4 


4a^ 


/ 1 


V29/2a2 


531 


(E7) 
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Here we used the fact that oc |a|^ and $o « 1 for |a| ^ 1. Thus we have still to calculate the asymptotics of iSap 
and iSsp which can be done by splitting the integrals, 


and 


o-3s 


1 /■ e 


1 


/ oo gS 2 r°° 

ds ——— exp[—4a^ cosh2sl-/ (ise“® exp[—4a^ cosh 2sl 

-oo cosh 2s ^ TT{8a-^)J_^ ^ 


27r(8a2) cosh 2 s 

|a|<i 1 Ki/2{4:a‘^ 

- 29/2^2 7r(8a2) 

|a|<l 1 1 


29/2^2 


+ 0{a-^) 


(E8) 
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-oo cosh 2s J ^(8a2)7_^ J 
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1 -ff3/2(4a2) 

■29/2„2+ ^( 3 ^ 2 ) 
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29/2^2 29/20fa3 213/20fa5 
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(E9) 


which yields the limit (104). 
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